Anderson localization of one-dimensional hybrid particles 
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We solve the Anderson localization problem on a two-leg ladder by the Fokker-Planck equation 
approach. The solution is exact in the weak disorder limit at a fixed inter-chain coupling. The 
study is motivated by progress in investigating the hybrid particles such as cavity polaritons. This 
application corresponds to parametrically different intra-chain hopping integrals (a "fast" chain 
coupled to a "slow" chain) . We show that the canonical Dorokhov-Mello-Pereyra- Kumar (DMPK) 
equation is insufficient for this problem. Indeed, the angular variables describing the eigenvectors of 
the transmission matrix enter into an extended DMPK equation in a non-trivial way, being entangled 
with the two transmission eigenvalues. This extended DMPK equation is solved analytically and the 
two Lyapunov exponents are obtained as functions of the parameters of the disordered ladder. The 
main result of the paper is that near the resonance energy, where the dispersion curves of the two 
decoupled and disorder-free chains intersect, the localization properties of the ladder are dominated 
by those of the slow chain. Away from the resonance they are dominated by the fast chain: a local 
excitation on the slow chain may travel a distance of the order of the localization length of the fast 
chain. 

PACS numbers: 72.15.Rn, 71.36.+C, 72.70.+m, 73.23.-b 



I. INTRODUCTION 

Despite more than half a century of history, Anderson 
localization 1 is still a very active field whose influence 
spreads throughout all of physics, from condensed mat- 
ter to wave propagation and imaging. 2 A special field 
where most of rigorous results on Anderson localization 
have been obtained, are one-dimensional and quasi-one- 
dimcnsional systems with uncorrelated disorder. Most of 
the efforts in this direction were made to obtain the statis- 
tics of localized wave functions in strictly one-dimensional 
continuous systems 4 ' 5 or tight-binding chains (see the re- 
cent work Ref. 22 and references therein). Alternatively, 
the limit of thick multi-channel N ^ 1 wires has been 
studied by the nonlinear super-symmetric sigma-modcl. 6 

A transfer matrix approach which allows one to con- 
sider any number of channels N was suggested by Gert- 
senshtein and Vasil'ev in the field of random waveguides. 7 
This approach has been applied to the problem of An- 
derson localization by Dorokhov 8 and later on by Mcllo, 
Pereyra and Kumar (DMPK). 12 It is similar in spirit to 
the derivation of the Fokker-Planck equation (the dif- 
fusion equation) from the Langevin equation of motion 
for a Brownian particle. However, in the present case 
an elementary step of dynamics in time is replaced by 
the scattering off an "elementary slice" of the Af-channel 
wire. As a result, a kind of Fokker-Planck equation 
arises which describes diffusion in the space of param- 
eters of the scattering matrix M, in which the role of 
time is played by the co-ordinate along the quasi-one- 
dimensional system. Usually the scattering matrix M is 
decomposed in a multiplicative way by the Bargmann's 
parametrization 12 which separates the "angle variables" 
of the [/(iV)-rotation matrices and the N eigenvalues 
T p =i,...,N of the transmission matrix. If the probability 



distribution of the scattering matrix is assumed invariant 
under rotation of the local basis (isotropy assumption), 
the canonical DMPK equation 8,11,12 may be obtained, 
which has the form of a Fokker-Planck equation in the 
space of N transmission eigenvalues. This equation was 
solved in Ref. 13 for an arbitrary number N of transmis- 
sion channels. 

The isotropy condition is not automatically fulfilled. It 
is believed that the isotropy condition is valid for a large 
number N ^> 1 of well coupled chains where the "elemen- 
tary slice" is a macroscopic object and the "local maxi- 
mum entropy ansatz" applies. 12 It is valid at weak disor- 
der in a strictly one-dimensional chain in the continuum 
limit a — > 0, or for a one-dimensional chain with finite 
lattice constant a outside the center-of-band anomaly. In 
this case the distribution of the only angular variable de- 
scribing a U(l) rotation, the scattering phase, is indeed 
flat. 22 

However, the case of few (N > 1) coupled chains is 
much more complicated. As was pointed out originally by 
Dorokhov 8 , and later on by Tartakovski, 10 in this case the 
angular and radial variables, are entangled in the Fokker- 
Planck equation. These are the variables determining the 
eigenvectors and eigenvalues of the transmission matrix, 
respectively We will refer to this generic Fokker-Planck 
equation as the extended DMPK equation in order to 
distinguish it from the canonical DMPK equation which 
contains only the radial part of the Laplace-Beltrami op- 
erator. The minimal model where such an entanglement 
is unavoidable, is the two-leg model of N = 2 coupled 
disordered chains. 

Yet this case is important not only as a minimal system 
where the canonical DMPK equation breaks down. It is 
relevant for the Anderson localization of linearly mixed 
hybrid particles such as polaritons. 33 Polaritons are the 
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result of coherent mixing of the electromagnetic field in 
a medium (photons in a waveguide for example) and ex- 
citations of matter (excitons) . In the absence of disorder 
photons have a much larger group velocity than excitons, 
and thus one subsystem is fast while the other one is 
slow. As a specific example, quasi-one-dimensional res- 
onators were recently fabricated by confining electromag- 
netic fields inside a semiconductor rod 14 or to a sequence 
of quantum wells. 15 In such resonators the dispersion of 
transverse-quantized photons is quadratic in the small 
momentum, with an effective mass as small as 1CT 4 of 
the effective mass of the Wannier-Mott exciton which is 
of the order of the mass of a free electron. 

Disorder is unavoidable in such systems due to the im- 
perfections of the resonator boundary and impurities. In 
many cases one can consider only one mode of transverse 
quantization for both the photon and the exciton. Thus a 
model of two dispersive modes (particles) with parametri- 
cally different transport properties arise. Due to the large 
dipole moment of the exciton these particles are mixed, 
resulting in avoided mode crossing. On top of that, dis- 
order acts on both of them, whereby its effect on the two 
channels can be rather different 16 . It is easy to see 16 that 
this system maps one-to-one onto a single particle model 
of two coupled chains in the presence of disorder. Ref. 16 
solved the coupled Dyson equations for the Green's func- 
tions of exciton and cavity phonon numerically, focussing 
on the so-called "motional narrowing" in the reflectivity 
spectra of normal incidence 18 . However, the issue of lo- 
calization of cavity polaritons was not raised. The latter 
was addressed in Ref. 17, which analyzed the scattering 
of electromagnetic waves in a disordered quantum-well 
structure supporting excitons. The random susceptibil- 
ity of excitons in each quantum well was shown to in- 
duce disorder for the light propagation, and the Dyson 
equation for the Green's function of the electromagnetic 
wave was then solved by the self-consistent theory of lo- 
calization. The author reached the conclusion that the 
localization length of light with frequencies within the 
polariton spectrum is substantially decreased due to en- 
hanced backscattering of light near the excitonic reso- 
nance. This is in qualitative agreement with our exact 
and more general study of the coupled disordered two- 
leg problem. The latter also finds natural applications 
in nanostructures and electronic propagation in hetero- 
geneous biological polymers, such as DNA molecules 19 . 

The main question we are asking in the present paper 
is: What happens to the localization properties when a 
fast chain is coupled to a slow one? Will the fast chain 
dominate the localization of the hybrid particle (e.g. a 
polariton) or the slow one? In other words: will the small- 
est Lyapunov exponent of the two-leg system (the inverse 
localization length) be similar to the one of the isolated 
fast chain, or rather to the one of the isolated slow chain? 
Can the presence of the "more strongly quantum" compo- 
nent (photon) help the "more classical" component (exci- 
ton) to get out of the swamp of localization? This latter 
question can be asked in many different physical situa- 
tions. It has been referred to as the "Munchhausen ef- 
fect" in Ref. 20, to describe the following effect predicted 



for a dc SQUID (superconducting quantum interference 
device) with two biased Josephson junctions, one with 
small plasma frequency (large mass), the other one with 
large plasma frequency (small mass): The junction with 
small mass can actually drag the "slower" junction (larger 
mass) out of its metastable state. 

Here we give an answer in the specific situation of a 
single hybrid particle. More interesting situations may 
arise when interacting and non-equilibrium polaritons are 
approaching Bose-condensation. 14,15 ' 36,37 ' 39 

A further question of more general interest can be ad- 
dressed by the same model problem. Namely, consider 
two or more coupled channels with similar propagation 
speed (i.e., inverse effective mass), but different disorder 
level: Which channel will dominate the localization, the 
cleaner or the more disordered one? This type of question 
arises not only in these hybrid single particle problems, 
but is an important element in the analysis of many par- 
ticle problems, where few and many particle excitations 
have various channels of propagation (e.g., all particles 
moving together, or moving in subgroups of fewer parti- 
cles). It is an important, but scarcely understood ques- 
tion, what determines the character of the propagation of 
such excitations when many parallel, but coupled chan- 
nels with different transport characteristics exist. Intu- 
itively one expects the fastest, and least disordered chan- 
nel to dominate the derealization. 

However, our analytical solution of the hybrid two-leg 
chain shows that in the one-dimensional case, this intu- 
ition is not always correct. Instead we find that, when 
the channels are strongly mixing with each other, it is the 
largest rate of back scattering, i.e., the more disordered 
chain, which dominates the physics. This may be seen 
as one of the many manifestations of the fact that in one 
dimension the localization length is essentially set by the 
mean free path. Our solution of the two chain problem 
furnishes a useful benchmark for approximate solutions 
in more complex and interacting situations. However, we 
caution that the phenomenology may be quite different 
in higher dimensions. We will discuss this further in the 
conclusion. 

The answer to the above questions will be obtained 
analytically from the exact solution of the two-leg (two- 
chain) Anderson localization model. This solution rep- 
resents a major technical advance, because for the first 
time a model, which leads to an extended DMPK equa- 
tion with non-separable angular and radial variables, is 
exactly solved. Without going into details our results are 
the following: 

(i) The answer depends qualitatively on whether or not 
the system is close to the resonance energy Er 1 which is 
defined as the energy where the dispersion curves of the 
two corresponding decoupled disorder-free chains inter- 
sect (see Fig. 1). 

(ii) Near the resonance the presence of the fast leg does 
not help to substantially delocalize the slow component 
(see Fig. 8). The localization length of a hybrid particle 
is at most by a factor of w 3 larger than the one of the 
slow particle (see Eqs. (95) and (96a)), being paramet- 
rically smaller than that of the fast particle. Thus the 
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slow particle dominates the localization properties of the 
hybrid particle near the resonance energy Er. 

(iii) A particular case where the resonance happens at 
all energies is the case of two coupled identical chains 
subject to different disorder (see Fig. 7). In this case the 
dominance of the more disordered chain extends to all 
energies thus pushing the localization length of the lad- 
der sharply down compared to that of the less disordered 
isolated chain. 

(iv) Away from the resonance the wavefunctions stay 
either mostly on the slow leg, being strongly localized. 
Or they have their main weight on the fast leg, and hy- 
bridize here and there with the slow leg (see Fig. 1 and 
Fig. 10). It is this second type of wavefunctions which 
helps excitations on the slow leg to delocalize due to the 
presence of the faster leg, even though this happens with 
small probability far from the resonance. 

(v) A very peculiar behavior occurs near the band-edges 
of the slow particle, where the system switches from two 
to one propagating channels. Just below the band-edge 
the localization length of the hybrid particle decreases 
dramatically being driven down by the localization length 
of the slow chain that vanishes at the band-edge (ne- 
glecting the Lifshitz tails). Above the band-edge the lo- 
calization length of a hybrid particle sharply recovers, 
approaching the value typical for the one-chain problem. 
Thus near the band-edge the localization length of the 
two-leg system has a sharp minimum, which is well re- 
produced by direct numerical simulations (see Fig. 9). 

The paper is organized as follows. In Section II the 
problem is formulated and the main definitions are given. 
In Section III the extended DMPK equation is derived. In 
Section IV the exact solution for the localization lengths 
is given and the main limiting cases are discussed. In 
Section V numerical results concerning the wave func- 
tions in each leg are presented. In Sec. VI a problem of 
one propagating channel and one evanescent channel is 
considered. The application of the theory to hybrid par- 
ticles such as polaritons, as well as considerations about 
higher dimensions, are discussed, in the Conclusion. 

II. TWO-LEG ANDERSON MODEL AND 
TRANSFER MATRIX FOR "ELEMENTARY" 
SLICE 

A. The model 

The Anderson model on a two-leg ladder is determined 
by the tight-binding Hamiltonian 

H = ^ ^ {e xv c\. v c xv - t v {c xv c x+lv + h.c)) 

i/=l,2 x 

X X 

where x e Z is the co-ordinate along the ladder, and 
v e {1,2} is the index labelling the two legs. In this 
model the on-site energies e xv are independently dis- 
tributed Gaussian random variables with zero mean, and 



t v is the hopping strength between nearest-neighbor sites 
on the v-th. leg. In general, the two legs will be subject 
to different random potentials, characterized by the two 
variances: 

°l = (2) 

We also consider different hopping strengths, for which 
we assume 

ti > *2- (3) 

The transverse hopping strength between the legs is t. 
Finally, it is natural to consider a homogeneous potential 
Se (i.e., a detuning) on leg 2. 

The Hamiltonian (1) is a generic model describing 
two coupled, uniformly disordered chains. Moreover, the 
model can also be adopted as an effective model to de- 
scribe non-interacting excitations with two linearly mix- 
ing channels of propagation in the presence of disorder. 
An important example is polaritons, the two channels 
correspond to the photon mode and the exciton mode, 
respectively. 

The model (1) has been studied analytically previously 
in the literature, focusing on the special case t\ = ti and 
o\ = a\. The continuous limit was solved long ago by 
Dorokhov. 8 The tight-binding model was considered later 
on by Kasner and Wcllcr. 9 Their results will be reference 
points for our more general study in the present work. 

The Schrodinger equation of the Hamiltonian (1) at a 
given energy E has the form 

*(a; - 1) + *(a: + 1) = (h(E) + e x )*(x), (4) 

where H>(x) is a single particle wavefunction with two 
components, representing the amplitudes on the leg 1 and 
2, 

HE) = (~f ~i e ) , (5) 

and 

6 ,=diag(^,^). (6) 

The terms h(E) and e x can be considered as the disorder- 
free and disordered part of the local Hamiltonian at the 
co-ordinate x. Notice that the disordered part (6) is ex- 
pressed as an effective disorder on the two legs, i.e., it 
is measured in units of the hopping strengths. In the 
analytical part of the present work, following the Fokker- 
Planck approach, we solve the problem exactly in the case 
of small disorder, He^H <C 1. 

1. Disorder free part 

The disorder-free ladder can easily be solved by di- 
agonalizing h(E) in Eq. (5). Thereby, the Schrodinger 
equation transforms into 

*(x-l) + *(x + l) = (h + e x )*{x), (7) 
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where 

h = diag(Ai, A a ), (8) 
and the "rotated" disorder potential is given by: 



■x+ 



e x _ cos 7 
_ sin 7 



e. T _ sin 7 
e x+ - e x - cos 7 



(9) 



Both depend implicitly on E via X T (E) and "f{E). In 
Eqs. (8) and (9) the following definitions are used: 
(i) In the disorder- free part (8), 
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FIG. 1: Two situations of clean energy dispersions. The 
dashed and solid curves correspond to decoupled and cou- 
pled chains. The decoupled dispersion curves intersect at the 
resonance energy Er- (a) no gap: E^ < E% ■ There are two 
propagating channels at a given energy for E± < E < E^- 
(b) gapped: E^ > E% ■ Apart from a forbidden band, only 
one propagating channel exists. 



where r G {1,2} is the channel or band index. As we will 
see in Eq. (15), r = 1 labels the conduction band, and 
t = 2 the valence band of the pure ladder. 38 
(ii) In the disordered part (9), 



Cxi 



£xl _j_ 6x2 



(11) 



are the symmetric and anti-symmetric combination of the 
disorder on the two legs. The "mixing angle" 7 = j(E) 
is defined through 



tan7(.E) = 



y/2ty/t\ + t 2 
(t 1 -t 2 )(E-E R y 



with a resonance pole at 



E R = Se — . 

h - h 



(12) 



(13) 



The value of 7 is chosen as: 7 € [0,7r/2] if E > Er] 
7 G [7r/2,vr] HE<E R . 

The pure system can be solved easily. In the absence 
of disorder the eigenfunctions ^(x) at energy E are com- 
posed of plane waves with momenta k T satisfying 



2 cos k T = X T . 



(14) 



±k T are degenerate solutions of Eq. (14), which is due 
to the space-inversion symmetry along the longitudinal 
direction of the pure ladder. Eq. (10) and (14) determine 
the energy dispersions of the conduction band and the 
valance band, 



E T (k) 



(ti +t 2 )cosk + 



Se 
~2 



(15) 



(ti - t 2 ) cosk + 



y 



t 2 . 



Generally, if ti ^ t 2 , the two decoupled bands (i.e., t = 
in Eq. (15)) cross at the energy Er [cf. Fig. 1], if \5e\ < 
2(<i — t 2 ). When the energy E is close to the resonance 



energy Er, the two legs mix with almost equal weights, 
even if we turn on a very small inter-chain coupling t. In 
the particular case of equal chain hoppings t\ = t 2 and 
no detuning Se = 0, there is a resonance at all energies 
since the two decoupled bands coincide. 

The top (+) and bottom (— ) edges of the r-band are 



E$=± (ti + ti) + 




<5e 

-t t ± 



(16) 



t 2 . 



According to Eq. (15), there are two cases of energy 
dispersions, which may arise depending on the choice of 
parameters: 

(i) In the case of E^ < E~^ [see Fig. 1(a)], there is 
no gap between the two bands. This is the case if the 
detuning Se and the interchain coupling t are both not 
too large. More precisely, one needs \Se\ < 2{t\ + 1 2 ) and 
t < t r , where 



t r . = 



4 (ti + t 2 f - Se 2 



tl+to 



(17) 



In the energy interval E^ < E < E^, we have two 
propagating channels; otherwise, at most one propagat- 
ing channel exists. 

(ii) In the opposite case, E± > E^ [see Fig. 1(b)], 
there is a gap between the two bands. We therefore have 
at most one propagating channel at any energy. 

Moreover, if k T is the wavevector of a propagating chan- 
nel, we call k T G (— 7r, tt], and k T > and k T < the 
right- and left-moving branch, resp. From Eq. (14) we 
also define a rapidity for each propagating channel as 



<9A T 
dk T 



(18) 
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2. Disordered part 

The impurity matrix (9) contains two ingredients which 
determine the localization properties of the model. One is 
e x ± [see Eq. (11)], which are the equally weighted (either 
symmetric or anti-symmetric) combinations of effective 
disorder on the two legs. The other is the mixing angle 7 
[see Eq. (12)], which describes the effective coupling be- 
tween the two legs. We refer to 7 as the bare mixing angle 
because it will be renormalizcd by disorder. The renor- 
malized mixing angle 7 [see Eq. (103)] will be discussed 
in Sec. IV. Being functions of these two quantities, the 
diagonal elements of e x are local random potentials ap- 
plied on the two channels r = 1,2, and the off-diagonal 
elements describe the random hopping between them. 

We analyze the model qualitatively in terms of effec- 
tive disorder and bare mixing angle before carrying out 
the detailed calculation. As already discussed above, ei- 
ther one or two propagating channels are permitted at a 
given energy. This leads to two distinct mechanisms of 
localization in the bulk of the energy band: 

(i) Two-channel regime. In this case, the physics is 
dominated by the mixing angle 7. If 7 ~ or n, the mix- 
ing of the two channels is weak: The magnitudes of off- 
diagonal elements of matrix (9) are much smaller than 
the magnitudes of the diagonal elements. This means 
that the two legs are weakly entangled, and the transverse 
hopping t can be treated as a perturbation. A perturba- 
tive study of wavefunctions in this regime is presented in 
Sec. VII. However, if 7 ~ it/2, the magnitudes of the 
off-diagonal elements are of the same order as the diago- 
nal elements. This implies that the two legs are strongly 
entangled. The localization properties are controlled by 
the leg with strong disorder, because in Eq. (11) it always 
dominates over the weaker disorder on the other leg. 

(ii) One- channel regime. The single-channel case has 
been solved by Bcrezinskii 4 and Mel'nikov 5 in the case of 
a single chain. The results they obtained can be applied 
in our problem by substituting the variance of disorder 
and rapidity of the propagating channel with the cor- 
responding quantities. However, we have to emphasize 
here that even if only one channel exists, coupling effects 
still present, since both the effective disorder and the ra- 
pidity in the remaining channel depend on the transport 
properties of both legs. In the one-channel regime the 
second channel is still present, but supports only evanes- 
cent modes. We show in Sec. VI that the effect of the 
evanescent channel on the propagating one is subleading 
when disorder is weak. 



B. Transfer matrix approach 

The Fokker-Planck approach and its related notations, 
such as the transfer matrix, the S-matrix etc. are intro- 
duced in detail in Refs. 11 and 12. We only outline the 
methodology here. The Fokker-Planck approach to onc- 
or quasi-one dimensional systems with static disorder at 
zero temperature is based on studying the statistical dis- 



tribution of random transfer matrices for a system of fi- 
nite length. An ensemble of such transfer matrices is con- 
structed by imposing appropriate symmetry constraints. 
In the present model there are two underlying symme- 
tries: time-reversal invariance and current conservation, 
which dramatically reduce the number of free parame- 
ters of transfer matrices. After a proper parametrization, 
the probability distribution function of these parameters 
completely describes the ensemble of transfer matrices, 
and therefore totally determines the statistical distribu- 
tions of many macroscopic quantities of the system, such 
as the conductance, etc. In order to obtain the probabil- 
ity distribution function of the free parameters, a stochas- 
tic evolution-like procedure is introduced by computing 
the variation of the probability distribution function of 
these parameters in a "bulk" system as an extra impu- 
rity "slice" is patched on one of its terminals, under the 
assumption that the patched slice is statistically inde- 
pendent of the bulk. Thereby, we construct a Markovian 
process for the probability distribution function. This is 
described by a kind of Fokker-Planck equation in the pa- 
rameter space of the transfer matrix with the length of the 
system as the time variable. Essentially, this procedure 
is analogous to deriving the diffusion equation from the 
Langevin equation for a Brownian particle. In practice, 
taking the length to infinity, we can analytically extract 
asymptotic properties of the model, such as localization 
lengths, etc., from the fixed point solution of the Fokker- 
Planck equation. 

As discussed above, the only microscopic quantity 
needed in order to write down the Fokker-Planck equation 
of our model Hamiltonian (1) is the transfer matrix of an 
"elementary slice" at any co-ordinate x. The Schrodingcr 
equation (7) can be represented in the following "transfer- 
matrix" form: 

$(x + l) =m x $(x), (19) 

where the 4-component wave function <f>(a;) and the 4x4 
transfer matrix m x is explicitly shown in the 2x2 "site- 
ancestor site" form: 

•«-(#(%)■ a-Cr--. 1 )- < 2 °> 

with ^/(x) and h, e x being the 2-component vector and 
2x2 matrices in the space of channels as defined in 
Eq. (7). 

The transfer matrix is manifestly real (which re- 
flects the time-reversal symmetry) and symplectic (which 
reflects the current conservation): 

m T x Jrh x = J, (21) 

where J is the standard skew-symmetric matrix: 

'-(-°iS)- <»> 

Note, however, that the transfer matrix rh x is not a 
convenient representation to construct a Fokker-Planck 
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equation. The reason is simple: because it is not diag- 
onal without impurities, the perturbative treatment of 
impurities is hard to perform. The proper transfer ma- 
trix m x is a certain rotation, which does not mix the two 
channels of the matrix m. x , but transforms to a more con- 
venient basis within each 2-dimensional channel subspace 
[see Appendix A]. The latter corresponds to the basis of 
solutions to the disorder-free Schrodinger equation ip T (x) 
(r = 1, 2 labeling the channels) which conserves the cur- 
rent along the ladder: 

j x = -i [i>*(x)il) T (x + 1) - h.c] = const = ±1. (23) 

For propagating modes with real wave vectors k T these 
are the right- and left- moving states 



which obey the conditions: 



(24) 



(25) 



For the evanescent modes with imaginary k = in the 
corresponding current-conserving states obeying Eq. (25) 
can be defined, too: 



± exp[q=«7r/4 — k t x] + exp[±i7r/4 + k t x] 

\/4 sinh k t 



(26) 



In this new basis of current-conserving states, the transfer 
matrix takes the form (see Appendix A): 



m x = 1 + 6m x , Sm x 



lCX.~E x OL a 



l x *z x u. x lOL x E x OL x 

ia x e x a x ia x e x a* x 



where the matrix ol x is diagonal in channel space, 



diag^+Oc),^^)) 



ch • 



(27) 



(28) 



Note that it is expressed in terms of the two components 
of the current conserving states Eqs. (24, 26) correspond- 
ing to the first and the second channel. 

In Eq. (27), the unit matrix 1 is the pure part of 
m x , which keeps the two incident plane waves invariant, 
and Sm x describes the impurities, which break the mo- 
mentum conservation and induce intra- and inter-channel 
scattering. The physical meaning of can be under- 
stood from the scattering processes described below. If 
there is only one right-moving plane wave in the 1-channel 
on the left-hand side (l.h.s.) of the slice, which is repre- 
sented by a four dimensional column vector with the first 
component one and the others zero, we can detect four 
components on the right-hand side (r.h.s.) of the slice, 
including the evanescent modes. In the case of two prop- 
agating channels these four components are right- and 
left-moving plane waves in the 1- and 2-channel, whose 
magnitudes and phase-shifts form the first row of m^. 
The other rows can be understood in the same manner. 
In short, the 11-, 12-, 21- and 22- block of ^m^ rep- 
resent respectively the right-moving forward-scattering, 
right-moving backward-scattering, left-moving backward- 
scattering and left-moving forward-scattering on the slice. 



In each block, the diagonal elements represent intra- 
channel scattering and the off-diagonal elements repre- 
sent inter-channel scattering. 

It is important that m x , Eq. (27), fulfills the same 
constraints regardless of the propagating or evanescent 
character of the modes [sec Appendix A] : 



m* = S 1 m x S 1 , mlS 3 m x = S 3 



(29) 



where S\ and S 3 are the four dimensional generalization 
of the first and third Pauli matrix with zero and unit en- 
tries replaced by 2 x 2 zero and unit matrices in the chan- 
nels space. The first condition follows from m* = m x , 
while the second condition is a consequence of the sym- 
plecticity Eq. (21). Thus these conditions are a direct 
consequence of the fact that m x belongs to the symplec- 
tic group Sp(4,R). As is obvious from the choice of the 
basis (23-26), their physical meaning is the time-reversal 
symmetry and the current conservation. 

The representation Eq. (27) of the transfer matrix of an 
"elementary slice" renders both the physical interpreta- 
tion and the symmetry constraints very transparent, and 
it will be seen to be a convenient starting point to con- 
struct the Fokker-Planck equation. On the other hand, 
since m,,, [see Eq. (20)] is real and has a relatively simple 
form, it is more suitable for numerical calculations. 



III. FOKKER-PLANCK EQUATION FOR THE 
DISTRIBUTION FUNCTION OF PARAMETERS 

A. Parametrization of transfer matrices 

Once the "building block" (27) is worked out, we can 
construct the Fokker-Planck equation by the blueprint of 
the Fokker-Planck approach. 5 > 8 > n > 12 The transfer matrix 
of a disordered sample with length L is 



M(L) = Y\_ mi = m L m L _i • mi, 



(30) 



x=l 



which is a 4 x 4 complex random matrix. It is easy to ver- 
ify that M.(L) also satisfies the time reversal invariance 
and current conservation conditions (29). It has been 
proved in Ref. 11 that all the 4x4 matrices satisfying 
Eq. (29) form a group which is identified with the sym- 
plectic group Sp(4,R). By the Bargmann's parametriza- 
tion of Sp(4,M), 12 one can represent M(L) as 



M 




u 

u* 



(31) 



where u and u are elements of the unitary group U(2), 
and statistically independent from each other, and 



F = diag(iq,F 2 ), 



(32) 



with F e e [l,oo) and g G {1, 2}. Because U(2) 21 has four 
real parameters, the group Sp(4, R) has ten real parame- 
ters. Furthermore, it is convenient to parametrize a U(2) 
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matrix by three Euler angles and a total phase angle, i.e., 

u(<f>, if, 9, ij>) = e-^e-^^e-^e-^*^ (33) 

in which a 2 and 173 are the second and third Pauli matrix, 
and the four angles take their values in the range <p, ip G 
[0, 2tt), 9 G [0, 7r) and tjj G [0,47r). In matrix form in the 
channels space, u can be written as 



u = e 2 



it /cos^e-5^+^) -sin I e~ 2 (¥>-</>) 



.5(V-^) 



, (34) 



ch 



which is convenient for the perturbative calculation be- 
low. The U(2) matrix u can be parametrized indepen- 
dently in the same form as Eq. (34). 

The probability distribution function of these ten real 
parameters determines completely the transfer matrix en- 
semble of the ladder described by the Hamiltonian (1). 
The goal of the Fokkcr-Planck approach is to obtain the 
Fokker-Planck equation satisfied by this probability dis- 
tribution function, in which the role of time is played by 
the length L. 

From Eq. (31) we obtain the transmission matrix 



{Ml 



u 



-1/2 



u. 



(35) 



by a simple relation between the transfer matrix and its 
corresponding S-matrix. 11,12,24 Due to the unitarity of u, 
the transmission co-efficients of the two channels are the 
two eigenvalues of the Hcrmitian matrix 



ttt 



u 



u 



t 



which are 



F + l 



(36) 



(37) 



where g G {1,2} is the index of the two-dimensional 
eigenspace of the matrix T. Now the physical meaning of 
the parametrization (31) becomes clear. The F g 's are re- 
lated to the two transmission co-efficients by the simple 
form Eq. (37). The matrix u diagonalizing the matrix T, 
contains the two eigenvectors of T, describing the polar- 
ization of the plane wave eigenmodes incident from the 
l.h.s. of the sample. For instance, if 9 = (u is a di- 
agonal matrix of redundant phases) , the two channels do 
not mix, and the incident waves are fully polarized in the 
basis of channels. On the other hand, if 6 = ir, the two 
channels are equally mixed, and the incident waves are 
unpolarized. In analogy to spherical co-ordinates, we will 
refer to the F e 's as the radial variables, while the angles 
in u or u are called angular variables. 

In principle, using the "building block" (27) and the 
parametrization (31), we can solve the full problem by 
writing down a Fokker-Planck equation for the joint prob- 
ability distribution function of all the ten parameters of 
M. However, since we are merely interested in the trans- 
mission co-efficients which arc determined by the proba- 
bility distribution function of T, instead of manipulating 



M, we study 
R MM' 



u 0\/ F VF 2 - l\ /ut 
u* ) { VW^l F ) I u T 



(38) 



R is a Hermitian matrix and contains only six parame- 
ters: 



X(R) = (F lt F 2 ,0,^,M- 



(39) 



The probability distribution function of A, denoted by 
Pl (A) , determines the transmission properties of the sam- 
ple with length L. Pl(A) is defined by 



Pi(A) = < y(X-A(R(L))) ) 



(40) 



where the overline denotes the average over realizations 
of the random potentials in the sample. It is convenient 
to introduce the characteristic function of Pl(A) 



P L [p) = J d\< >l! - X P L [\) = r'/'-AiltlL!). 



(41) 



Our main goal in this paper is to calculate the two 
localization lengths, defined as the inverse Lyapunov ex- 
ponents of the transfer matrix (31): 



1 d 



e i(2) = _ i^2^ lnTm ^ (min))i ' 



(42) 



in which the subscripts "max" and "min" denote the 
larger and smaller of the two real values Ti^, and the 
averaging (-) L is earned out with the probability distri- 
bution Pt (A). Therefore, by definition 



£1 > 6- 



B. Physical interpretation of £1 and £2 



(43) 



It is worthwhile to visualize how two parametrically dif- 
ferent localization lengths (£1 3> £2) manifest themselves 
in transport properties. For instance, let us discuss the 
dimensionless conductance g = Ti + T 2 , a typical behav- 
ior of which is shown as a function of the sample length 
L in Fig. 2. If L <C £2, T e ~ 1, and g ~ 2 corresponds to 
a nearly perfect transmission. As L increases, T e and g 
decay exponentially. T 2 decays much faster than Ti since 
£1 ^ £2- As long as L < £1 the system still conducts 
well since g is still appreciable. For L ~ £1 it crosses over 
to an insulating regime. On the other hand, £2 marks 
the crossover length scale below which g(L) (black curve) 
decreases as fast as T 2 (blue dashed curve) until the con- 
ductance saturates to a plateau g w 1. For L > £ 2 , g 
decays with the slow rate £f , like Ti (red dotted curve). 
Therefore, the two parametrically different localization 
lengths can be identified by two distinct decay rates of g 
at small and large length scales. 
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FIG. 2: Schematic diagram for typical values of dimensionless 
conductance g as a function of the length L in the case £1 
£2- A crossover happens at L ~ £2- In the region L < £2, g 
decays as fast as Tmm (see the insert). Once L > £2, g decays 
as slowly as T max . The system is well conducting if L < £1 
and crosses to an insulating regime for L > £1. 



The statistics of transmission eigenvalues and lo- 
calization lengths of disordered multi-channel micro- 
waveguides have been visualized in experiments . 25 How- 
ever, only more or less isotropically disordered cases 
(identical hopping and disorder strength in each chan- 
nel) were realized, while a situation where £1 ^> £ 2 is 
hard to achieve in such systems (see Ref. 25 and refer- 
ences therein). In contrast such anisotropic situations 
are rather natural in exciton polariton systems. 

We will see in Sec. VII that the two localization lengths 
£1 and £2 also characterize the spatial variations of the 
eigenfunctions ^(x) on the two legs. 



C. Fokker-Planck equation for the distribution 
function of A(R) 

Having a disordered sample of length L, whose transfer 
matrix is M(L) and adding one more slice, we obtain the 
transfer matrix of the sample with length L + 1: 



M(L + 1) = m L+1 M(L). 



(44) 



Simultaneously, according to Eqs. (38) and (27), R(L) is 
updated to 

R(L + 1) = m L+ iR(L)m[ +1 = R(L) + SR, (45a) 



SR = (R5mjr +1 + h.c.) + 6m L+ iR<5m[ +1 . (45b) 
Accordingly, A(R(L)) is incremented by 

X(R(L + 1)) = X(R(L)) + SX. (46) 



According to Eqs. (41) and (46), we obtain the char- 
acteristic function of Pl + i(X): 



P L+1 (p) = e iP'A(R(L+l)) = e iff-X(R(L)) e ip-8X_ 



(47) 



We can expand e l ^' sx on the r.h.s. of Eq. (47) into a 
Taylor series e& sX = J2n i^' S *) n - Usin S E( l s - (45,46) 
standard perturbation theory yields an expansion of <5A 
in powers of the disorder potential as SX — J2 n >i SX^ n \ 
where SX^ is of n-th order in e. With this, the r.h.s. of 
Eq. (47) can be expanded in powers of the disorder po- 
tential. In principle, we can proceed with this expansion 
to arbitrarily high orders. Thereafter, the average over 
disorder on the slice L + 1 can be performed. Eqs. (45)- 
(47) fully define our problem. However, it is impossible 
to solve it analytically without further simplification. 

Progress can be made by considering the weak disor- 
der limit. In the two-channel regime, the weak disorder 
limit implies that both of the mean free paths are much 
larger than the lattice constant. As a first estimation, ap- 
plying the Born approximation to an "elementary slice" , 
the inverse mean free paths of the two propagating chan- 
nels can be expressed as certain linear combinations of 
the variances of the effective disorders on the two chains, 
denned as 



xl 



(48) 



for chain v. In the weak disorder limit where the smaller 
of the two localization lengths is much larger than the 
lattice constant, 



I > 1, 



(49) 



only the terms proportional to xl on t ne r.h.s of Eq. (47) 
have to be taken into account. Hence, we calculate SX 
perturbatively up to the second order [see Appendix B]. 
If L >• 1, as we always assume, Pl+i — Pl — 9lPl- 
Under these conditions, Eq. (47) leads to 



d L P L =ip-8X( 2 )e^ x ^{L)) 



V 



■ SXW) e v-MR(£)). 

Note that because the random potentials in different 
slices are uncorrelated, the terms can be averaged 

independently of e 1 ?' A ( R ( L )) . By the inverse of the Fourier 
transform defined in Eq. (41) we obtain the Fokker- 
Planck equation for Pl(X): 



d L P 



i=i 



(2) 



In Eq. (51) the averages are taken over the realizations 
of random potentials in the slice at L + 1. 

The Fokker-Planck equation (51) can be rewritten in 
the form of a continuity equation: 



d L P 



(52) 
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where the generalized current density Jj takes the form: 

6 

Ji = Vi (A)P - D H W 9 ^ P > ( 53 ) 

with 

v i (\) = 5)$ l) +dx j D ij (\), (54a) 



D ll (\) = -5\ { p5\f\ (54b) 

Vi(A) and Dij(X) are a generalized stream velocity and a 
generalized diffusion tensor, respectively. 

In order to solve Eq. (51), we have to add the initial 
condition, namely, Pq(X). Usually Po(X) is chosen as the 
probability distribution function in the ballistic limit, 13 

P (A) = S(F 1 - 1)S(F 2 - l)6(6)6(<f>)6(ij)6(ip), (55) 

where b(x) is the Dirac delta function. However, as we 
will see later, a unique fixed point of Pl(X) exists in the 
limit L — >• oo, which does not depend upon the initial 
condition. Essentially, the existence of a fixed-point so- 
lution of Eq. (51) is protected by Anderson localization 
which prevents the system from chaos. 26 

D. Coarsegraining 

Let us analyze the r.h.s. ofEq. (51) qualitatively. From 
Eqs. (BIO) and (B12) in the Appendix, it is clear that the 

co-efficients SX^SX^ and bX^ are sums of terms carry- 
ing phase factors 1, e ±*( fc i- fe 2)£ j _ a nd so on These phase 
factors come from the disorder average of products of two 
elements of the matrices (27). Their phases correspond to 
the possible wave vector transfers of two scatterings from 
a slice, similarly as found in the Berezinskii technique 4 . 
They are thus linear combinations of two or four values 
of ±ki i2 : 

K osc = {±Afc, ±2Afc, (56) 
±2fc 1(2 ), ±(k 1 +k 2 ), ± [3fc 1(2) - fc 2 (i)] , 
±4fc 1(2) , ±2(k 1 + k 2 ), ± [3fci (2) + k 2{1) ] } , 

where 

Ak = fa-k 2 . (57) 

Terms with phase "0" do not oscillate. The largest spatial 
period of the oscillating terms is 

-^osc = max Sk^ 1 (58) 
<5fceK OB( . 

Under the condition that 

Lose < I, (59) 

a coarsegrained probability distribution function can be 
defined as the average of Pl(X) over Lose- From now on, 



we use the same symbol Pj, (A) to denote its coarsegrained 
counterpart, which satisfies Eq. (51), but neglecting the 
oscillating terms. 

Additionally, at special energies it may happen that an 
oscillation period becomes commensurate with the lattice 
spacing , 5k = n/n. An important example of this com- 
mensurability is the situation where Se = 0, 2(k\ + k 2 ) = 
2tt at E = 0. In this case the terms with the phase 
factor e ±2l ( k i+ k 2)L J,-, n0 ^ avera g e a nd give anomalous 
contributions to the non-oscillating co-efficients. This ef- 
fect leads to the so-called center-of-band anomaly in the 
eigenfunction statistics of the one-chain Anderson model 
(see Ref. 22 and references therein). While they are not 
included in our analytical study, the commensurability- 
induced anomalies can be seen clearly in the numerical 
results for localization lengths (cf. Figs. 6, 7, 10 and 11). 

The coarsegraining procedure leads to a significant sim- 
plification: the co-efficients on the r.h.s. of Eq. (51) do 
not depend on L, <fi arL d ip any longer, which renders 
the solution of Eq. (51) much easier. Its non-oscillating 
co-efficients are evaluated in Appendix B. We do not re- 
produce them explicitly here, since we further transform 
the Fokker Planck equation below. However it is worth- 
while pointing out a formal property of its co-cfficicnts. 
From Eqs. (B8), (BIO) and (B12), it is easy to see that 

the ingredients for evaluating bX^p bX^ and bX^ arc the 
disorder-averaged correlators between any two elements 
of matrices (27). During the calculation, three Born cross 
sections appear naturally, being covariances of the effec- 
tive disorder variables, 

^^(^cos^+^sin 4 !), (60a) 
y 2 = ^( Xl sin^+ X > S 4 !), (60b) 

^3 = i ^(x 1 +^)sin 2 ^cos 2 |, (60c) 

in which V\(2) corresponds to intra-channel scattering 
processes fci( 2 ) <h> — &i( 2 ), and U 3 corresponds to inter- 
channel scattering processes fc 1 ( 2 ) «-» ±k 2 ny Note that 
the effective disorder variances (48) enter into the three 
Born cross-sections, instead of the bare variances (2). We 
will see that the above three Born cross-sections com- 
pletely define the localization lengths and most phenom- 
ena can be understood based on them. 

We note that the coarsegraining, through Eq. (59), im- 
poses a crucial restriction on the applicability of the sim- 
plified Fokker-Planck equation. According to Eqs. (10) 
and (14), if t is small enough, at E = Er, 

|Afc|oci. (61) 

In this case, Eqs. (59) and (61) require that 

t > bE, (62) 
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where; 

SEoaT 1 , (63) 

is the characteristic disorder energy scale (essentially the 
level spacing in the localization volume). In other words, 
Eq. (62) imposes a "strong coupling" between the two 
legs, as compared with the disorder scale. However, from 
the point of view of the strength of disorder, Eq. (62) is a 
more restrictive condition than £ 2 3> 1 011 the smallncss 
of disorder. However, it is automatically fulfilled in the 
limit cr„ — > at fixed values of coupling constants t and 
t„. 

Eq. (62) restricts the region of applicability of the sim- 
plified equation (68) which we will derive below. Indeed, 
we will see that by simply taking the limit t — > in the 
solution of that equation one docs not recover the triv- 
ial result for the uncoupled chains. This is because the 
equation is derived under the condition that t is limited 
from below by Eq. (62). The "weak coupling" regime is 
studied numerically in Sec. V A and the cross-over to the 
limit of uncoupled chains is observed at a scale of t ~ l^ 1 
as expected. 

Since the definition of localization lengths (42) only 



involves the F e 'a, and since the co-efficients of Eq. (51) do 
not contain <j) and ip, we define the marginal probability 
distribution function 

W L {F l ,F 2 ,e,i>) = J d<j>dyP L {\). (64) 

Further we change variables to the set 

ff=(F u F 2 ,u,i;), (65) 

where 

m = cos6», we (-1,1]. (66) 

We thus have 

W L (ff} = J d0dpP L (Fi,F 2 ,0(u),lM,¥>). (67) 

Substituting Eq. (67) into (51), and replacing the differ- 
ential operators do — > —y/l — u 2 d u and dg — > — ud u + (l — 
u 2 )d 2 , we obtain the Fokkcr-Planck equation for Wl{vi)'- 



(68) 



i=l 



j>i=l 



The co-efficients q , are relatively simple functions of 
ff. They can be obtained from the averages of the matrix 
elements computed in App. B, and are given in App. C. 
However, only a small number of them will turn out to 
be relevant for the quantities of interest to us. 

One can see that in Eq. (68) the radial variables, 
F g , arc entangled with the angular variables u and ip. 
Thus Eq. (68) is more general than the canonical DMPK 
equation 8 ' 11,12 where only radial variables appear. To 
emphasize the difference we refer to Eq. (68) as the ex- 
tended DMPK equation. The derivation of Eq. (68) for 
the two-leg problem is our main technical achievement in 
the present paper. It allows us to obtain the evolution (as 
a function of L) of the expectation value of any quantity 
defined in //-space. 



IV. CALCULATING THE LOCALIZATION 
LENGTH 

It is well-known that in quasi-one dimensional settings 
single particles are always localized at any energy in ar- 
bitrarily weak (uncorrelated) disorder. 3 The localization 
length quantifies the localization tendency in real space. 
In this section we calculate the localization lengths for 
the present model. 

The analytic expression of In T max ( min ) in Eq. (42) can 



be written as 

lnT max(min) = e(AF)lnT 2(1) + 6(-AF) lnT 1(2) , (69) 

where Q(z) is the Heaviside step function and 

AF = Fi — F 2 . (70) 

Multiplying both sides of Eq. (68) by the r.h.s. of Eq. (69) 
and integrating over all the variables, we obtain from 
Eq. (42) 

C; 1 = lim ((A + D 2 ) - (-l)'(Di - D 2 )sgn(AF)) L , 

(71) 



L— >oo 



with 



Di = 



(Fi + iy 



+ 



dFiCu 



F + l Fi + 1 



, ie{l,2}, 



(72) 

in which p G {1,2}, sgn(2;) is the sign function and the 
co-efficients Ci,ca [see App. C] are 



(-1)<2- 



: r 6 , 



with 



AF 

Cll = {Ff - l)V u 

Yi{u) =V! +V 2 + W 3 + (-1Y2 (V 2 - V-l) u 
+ (V 1 + V 2 -W 3 )u 2 , 

r 6 (u) =vi + v 2 - (Vi + v 2 - w 3 )u 2 . 
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The formula (71) for the localization lengths can be fur- 
ther simplified in the limit L 3> 1. When L is large, the 
typical value of F min ( max ) is of the order of e L /^ 1 < 2 > , which 
is exponentially large. Therefore, 



^max -^min 1 > 



(73) 



as we assume £1 > £ 2 (see Eq. (43)). The hierarchy (73) 
largely simplifies the co-efficients of Eq. (68), which leads 
to 



Cl 



lim 

L^co Fi + 1 



c 2 



lim 

L^oo F 2 + 1 



lim 



L^oo (Fi + l) 2 



-2r 6 e(AF), 



-2r 6 e(-AF), 



= r, 



(74a) 



(74b) 



(74c) 



(74d) 



lim = 2T«. 

L^oo Fi + 1 

As a result, Eq. (71) reduces to 
C; 1 =Vi + V 2 + 2V 3 

+ (-iy( K l(T e ) + \v 1 -v 2 \(u)y ( ' r " 

where (•) = limL^ 00 (-)i, Vi, V 2 and V 3 are the Born 
cross sections defined in Eq. (60). The main simplifica- 
tion is that T 6 depends only on u, but not on the other 



parameters of the scattering matrix. Therefore, the lo- 
calization lengths are fully determined by the marginal 
probability distribution function of u defined by 

Wi(«) = J dFidFidipWLiff). (76) 



Integrating over F\ , F 2 and ip on both sides of Eq. 
we obtain the Fokker-Planck equation for wi(m): 

d L w = d u (c 33 d u w) + d u {c 3 w), 



(77) 



where c 3 and c 33 are derived in App. C. It has a fixed- 
point solution satisfying 

d u (c 33 d u w) + 8 U (c 3 w) =0. (78) 

In the large L limit the co-efficients are given by 

lim c 3 = (IV, -V 2 \- d u T 6 ) (1 - u 2 ), (79a) 



lim c 33 = (V 3 + r 6 )(l- U 2 ). (79b) 



From Eq. (79) one can see that in the limit (73), c 33 
and c 3 do not depend on F± , F 2 and ip any longer. There- 
fore, Eq. (78) is reduced to an ordinary differential equa- 
tion with respect to u. By considering the general con- 
straints on a probability distribution function, namely the 
non-negativity w(u) > and the normalization condition 
/ duw(u) = 1, the solution to Eq. (78) is unique, 



w(u) 



Qi cxp ( — arctan — ) 

wi(u) = — % U , AV<0, 



w 2 (w) = 



91 



92+* 
, 82-* 
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2sinh(^lnfa±i)(g!-„2)' 



(80) 



AV > 0, 



where 

Ql = 2\V 1 -V 2 \/\AV\, (81) 
92 = V(Vi + V 2 +W 3 )/\AV\, (82) 

and 

AV = V 1 +V 2 -W 3 . (83) 

Eq. (75) and (80) are our main analytical results. The 
localization lengths are expressed entirely in terms of the 
three Born cross sections V\, V 2 and V 3 . We recall that 
we made the assumptions of weak disorder, Eq. (49), and 
sufficiently strong coupling, Eq. (62). 

In Eq. (80), w 2 (u) is simply the analytical continuation 
of wi(«). To show this, we start form AV > side and 



drop the absolute value on AV. If AV crosses zero from 
above, namely AV — > —AV, qi changes continuously to 
—qi, and q 2 changes to one of the two branches ±z|g 2 | 
because of the square root. It can be easily verified that 

wi(u;-9i,±i|g 2 |) = w 2 (m;<7i, \q 2 \), (84) 

by the formula arctan z — i/21n[(l — iz)/(l + iz)] for a 
complex number z. 

Given the physical meaning of the parameter u, it is 
natural to interpret the analytical continuation as de- 
scribing the crossover between two regimes of the po- 
larization, as controlled by the relative strength of the 
effective disorders. If AV > (i.e., V\ + V 2 > 4V 3 ) the 
intra-channel scattering is stronger than the inter-channel 



12 



0.8 


Single- 
channel 






Single- 
\ channel 


0.6 










0.4 




Resonance 




0.2 








^Off- \ 


0.0 


/ resonance 






resonance \ 




-0.4 


-0.2 

E 


0.0 


0.2 



FIG. 3: The resonant and off-resonant regimes on the E — t 
plane for ti = 1, ti = 0.2, with isotropic disorder a\ — o\ 
and bias Se = —0.1. The black curve indicates the band edges 
E = and E = E% , beyond which only one channel exists. 
The blue curve marks the crossover line AV = 0. At small 
enough t, AV = can be linearized to t ~ k\E — Er\, cf. 
Eq. (85), which is plotted as the red dashed line. 



scattering, while AV < means the opposite. The two 
regimes can be distinguished quantitatively. According to 
Eq. (60), the co-efficients in the linear combination of the 
effective disorder parameters, x 2 ,, are determined by the 
bare "mixing angle" 7 and the rapidities, v T . Suppose the 
resonance energy En is approached while keeping t < t c 
[see Eq. (17)]. If E is in the vicinity of En, 7 ~ n/2, 
and < 0. Otherwise, if E is far enough from En, 
7 — > or 7r, and AV > 0. Therefore, there must be an 
energy interval around En, in which the physics is similar 
to that at resonance, 7 = 7r/2. Further away from En 
the physics is similar to the limiting cases 7 — or tt. We 
call AV^ < and AV^ > the resonant and off-resonant 
regimes, whose distinct behavior we will analyze below. 



A. Resonant and off-resonant regimes 

As shown in Fig. 3, for fixed t v and Se, AV = (blue 
curve) divides the E — t plane into two regions in the two- 



channel regime (below the black curve). Three important 
observations are in order: 

(i) At weak coupling t, more precisely, for (<f c , but 
still within the condition (62), the relation AV = for the 
border of the resonance region implies the linear relation 
(see the red dashed lines in Fig. 3) 



with 



t ~ n(ti,t 2 )\E - Ef 



n(t 1 ,t 2 ) 



(85) 



(86) 



The slope K{t\,t2) neither depend on cr 2 . nor on Se. 

(ii) If the coupling t is strong enough, the resonance 
energy interval shrinks to zero as t — > t c (the top edge of 
Fig. 3). This "re-entrance" behavior is due to the compe- 
tition between the strong coupling, which pulls 7 close to 
7r/2, and the band edge effect, which reduces the rapidity 
of one of the channels. We can illustrate this behavior by 
considering two limiting cases. If t is weak, its effect is of 
first order on 7, but of second order on the v T . Therefore, 
the coupling wins and the resonance energy interval fol- 
lows the linear relation (85). Alternatively, if the energy 
is in the vicinity of the band edges E = E± and E% , one 
of the rapidities tends to zero. As a consequence, V\ or V2 
is much larger than V3, which gives a large positive AV. 
Therefore, there is always some region around the band 
edges (black curves in Fig. 3), which is out of resonance. 
As the crossover line must match the two limits t — >• 
and t t c , it is necessarily re-entrant. 

(iii) In the case of a non-zero detuning energy Se the 
resonant energy interval is slightly asymmetric around 
E = E R . 



B. Fixed point distribution w(u = cos (9) 

Let us now discuss the distribution Eq. (80) in different 
regimes and some of its consequences. For this purpose, 
we plot in Fig. 4 some representative w(u) together with 
the expectation value and variance of u. We select various 
values of E across the resonant and off-resonant regime. 
Two types of behavior can be observed in the two regimes: 



(i) Near the resonance, u — cos 6 is distributed relatively 
uniformly in the interval (—1,1]. Its average value is 
much smaller than 1, but its variance is large of order 
O(l). However, the distribution is definitely not com- 
pletely uniform. Indeed, the limit of the distribution can 
be obtained form Eq. (80) in the weak coupling limit as 
t -> as 



which is manifestly non-uniform. A similar distribution 
was obtained by Dorokhov 8 in the case of two equivalent 
chains. We will discuss the difference to Eq. (87) later. 

(ii) Off resonance, the distribution function w(u) is 
strongly peaked at u = 1, and its fluctuations are strongly 
suppressed. 

At this point the difference between the canonical 
DMPK equation, 8,11 ' 12 which applies in the case N 3> 1, 
and the extended DMPK equation obtained here for the 
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FIG. 4: Marginal probability distribution for the angular variable u — cos((9), w(u), in different regimes. Here the inter-chain 
coupling t = 0.1 and the other parameters are the same as in Fig. 3. (a) Distributions w(it) for three energies across the 
resonance and off-resonance regimes, namely, E — Er rs —0.13 (i.e., resonance energy), E — 0.2 (off-resonance) and E = 0.0 
(crossover). On resonance the distribution is nearly uniform, while it is strongly non-uniform off-resonance (b) The expectation 
value and variance of u as functions of E. 



case N — 2, is clear. The isotropy assumption, which al- 
lows one to derive the canonical DMPK equation, states 
that the angular variable distribution w(u) should be uni- 
form, i.e., independent of u, in contrast to Eq. (87). In 
order to justify the canonical DMPK equation, we have 
to have a large number of equal chains. A sufficient con- 
dition for obtaining the canonical DMPK equation is that 
the probability distribution of the transfer matrices of an 
"elementary slice" is invariant under U(N) rotation. This 
situation may be achieved in thick wires. 11,12 However, in 
few-channel cases the localization lengths are larger, but 
still of the same order as the mean free path. There is 
no parametric window between them that permits the 
emergence of ?7(7V)-invariant ensembles of transfer ma- 
trices upon coarsegraining. 

The qualitative difference in the distribution function 
w(u) in the two regimes has important implications on the 
localization lengths. To calculate the localization lengths 
from Eq. (75), we need (r 6 ) and (u). Using Eq. (80) we 
obtain 



(r 6 



q 1 q2\AV\S 1 (q 1 ,q 2 ) 

2 sinh ( — arctan — 
V 12 12 , 

qig2|AV|gi(qi,g 2 ) 
2sinh(iln«±i) 

V 2 12 12 - 1 J 



4V 3 AV < 0, 



4TA AV > 0, 



and 



(u) 



qi Sa (giwa) 



2 sinh I 



qiS 2 (qi,q2) 



2 sinh 



(SL. 
\2l2 



In 



92 + 1 



AV < 0, 



AV > 0, 



(89) 



where Si(2) and <Si(2) are integrals defined by 



Si(qi,q2) = 

^2(91,92) = 



arctan (I/92) q± 

dz sec 2 ze~i2 z , 

— arctan (1/92) 
arctan (I/1J2) n 

dz tan z e q 2 , 

- arctan (1/52) 
1/92 

dz 

-1/92 V 1 - Z 
1/92 
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1 + Z \ 2q 2 



(90) 



dz- 



1/92 



1 - z 2 \ \-z 



91 

1 + Z \ 2 12 



C. Numerical analysis 

In order to confirm our analytical results for the local- 
ization lengths in Eq. (75) we calculated numerically the 
Lyapunov exponents of the products of transfer matrices 
in Eq. (30). An efficient numerical method, known as the 
reorthogonalization method, has been developed in the 
study of dynamical systems 31 and widely spread in the 
field of Anderson localization. 27 The forthcoming numer- 
ical results in Figs. 5-8, 10 and 11 are all obtained by this 
method. 

The usefulness of the reorthogonalization method is not 
restricted to numerical simulations. It also provides the 
basis for the perturbative analysis about the Lyapunov 
exponents in the weak disorder limit in Sec. VI B. 
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V. RESULTS FOR THE LOCALIZATION 
LENGTHS 

In order to reveal the effects of the transverse coupling 
t on the localization lengths, we define the two ratios 



r P =U6°\ Pe{l,2}, 



(91) 



Notice that in the resonant case, the Born cross- 
sections (60) are dominated by x!> which gives rise to 
the dramatic drop of the localization length of the fast 
leg: The slow-leg is dominating the backscattering rate 
and thus the localization length. 

From Eq. (98) we draw several important conclusions 
below. 



where the £p°^s are the localization lengths of the decou- 



pled legs, for which we may assume > For 
simplicity, we refer to the leg 1 and 2 as the fast- and the 
slow-leg, respectively. The bare localization lengths £p ' 
can easily be obtained from Eq. (75) by taking 7 = 0, 
w(u) = 5(u — 1) and t = 0, which yields 



,(0) 



t(0) 
So 



A.p 



(92) 



Eq. (92) coincides with the well-known single-chain 
result. 5 



A. E = and 8e — 0: Resonant regime 

Consider first the case 5e = 0, in which the resonance 
energy vanishes E R = 0. From Eqs. (12) and (18) it 
follows that the mixing angle is 7 = tt/2 once t 7^ 0, and 
the two rapidities v\ = v 2 — v are equal to each other: 



v 2 =4- 



(93) 



Consequently, the three Born cross-sections have the 
same value and are equal to: 



V 1 = V 2 = V 3 = V 



1 



16v 



2 {xl + xl)- (94) 



This gives q\ = and q 2 = V% according to Eqs. (81) 
and (83). Evaluating the integrals (90), we obtain the 
two localization lengths 



t P = 8C p v 2 /(xl+xl), 



where 



and 



Ci = 



3(tt - \/3) 



0.743, 



(95) 
(96a) 



1. Statistically identical chains 

For two coupled chains, which are statistically identi- 
cal, one has xi = x\-> an d we obtain 

-|y = n = 2d « 1.486, 
s"i 



AO) 

S2 



= r 2 = 2C 2 w 0.754. 



(99) 



We note that n is slightly larger than the value obtained 
by Dorokhov 8 , which is w/(w—l) w 1.467. The reason is 
that we have taken into account the forward-scattering in 
the "elementary slice" (27), which was neglected in the 
work by Dorokhov. Moreover the latter was restricted 
to t\ — t 2 - In Fig. 5 we compare our analytical predic- 
tion with Dorokhov's. Note that we take E — 0.1 in the 
numerical simulation in order to avoid the anomaly at 
E = 0, as mentioned in Sec. Ill D. The enhancement fac- 
tor n is essentially independent of the selected energy if 
t is weak enough. This is due to the fact that any energy 
is at resonance conditions for t\ = t 2 . 

The effect of forward scattering, which was included 
in our work, is clearly visible. It is confirmed by the 
numerical simulation at resonance conditions. However, 
the value r\ ~ 1.776 obtained by Kasner and Wellcr 9 
deviates significantly from our numerical and analytical 
results. 



2. Parametrically different chains 

It is interesting to analyze what happens if the bare 
localization lengths of the chains are parametrically dif- 



ferent £2°^ 



<C In the resonant regime, for W 2 , \E — 



E R \ <C t -C ti,t 2 , we obtain: 



(o) 



2.972 ^ 0) , 



(100a) 



C 2 = = w 0.377. 

7T + 3V3 



(96b) 



The corresponding decoupled values (t = 0) can easily be 
obtained from Eq. (92), 

Therefore, the ratios defined by Eq. (91) read 



(98) 



6^4C 2 ^ 0) « 1.507^ 



;(0) 



(100b) 



Eq. (100) is one of the central results in this paper: In 
the resonant regime, the localization length of the fast 
leg is dramatically dragged down by the slow-leg. In con- 
trast, the localization length of the slow leg is increased 
by the presence of the fast leg, but remains of the same 
order. As a result both localization lengths become of 
the order of that for the bare slow leg. This is illustrated 
for two different cases of coupled fast and slow legs in 
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FIG. 5: Ratio of coupled and uncoupled localization lengths, 
n = for statistically identical chains with ii = ti = 1 

and small disorder W. We consider small couplings t, for 
which any energy is at resonance conditions. In the numerical 
simulation we take E = 0.1 in order to avoid the anomaly at 
E = 0.0. However, according to Eq. (99), r\ is almost inde- 
pendent of E if t is weak enough. The continuum approxima- 
tion becomes exact for W 2 <C t <C 1, as is illustrated by the 
convergence of the numerical data to the analytical prediction 
(the agreement is already good for cW 2 < t with c « 0.25). 
For comparison we also plot Dorokhov's prediction 8 , which 
neglected forward scattering in the Fokker-Planck equation. 
The result of Kasner and Weller 9 (n w 1.776) is in clear con- 
tradiction with these numerics. 



Figs. 6 and 7. Fig. 6 shows the effect in the case of legs 
with equal disorder but different hopping strength, the 
resonance being at E — 0. In Fig. 7 the faster leg has 
the same hopping but weaker disorder. Here the legs are 
resonant at every energy below the band-edge . 

We note that there is no regime where both n > 1 and 
ri > 1, as this would contradict the equality r p /C p = 
v 2 , which follows from Eq. (98). At the band center and 
t — » one can achieve that both localization lengths do 
not decrease upon coupling the chains, r% = r% = 1. This 
happens when xi/x? = 4Ci — 1, which assures that the 
localization lengths do not change at coupling constants 
t < W 2 according to the discussion in Sec. VA3. 



3. Weak coupling limit 



Upon simply taking the t — limit, 

r P = 4C pX 2 p /(xl+xl), 



(101) 

one does not recover the decoupled values r p = 1. This 
should indeed be expected, as we have already discussed 
in Sec. IIIC. The reason traces back to condition (62) 
to obtain Eq. (68), namely that t be larger than the 
disorder energy scale SE oc W 2 . In order to verify the 
non-commutativity of t — > and W — > 0, we computed 
numerically the localization lengths by the transfer ma- 
trix approach, and obtained the values of r p down to very 



FIG. 6: Localization lengths for chains with different hopping 
strengths (ti = 1 and t% = 0.1), but equal disorder (W = 0.05) 
as a function of energy at detuning 8e = and intermedi- 
ate coupling t = 0.1. The solid curves are analytical results. 
Black curves correspond to uncoupled chains, red ones to the 
coupled chains. The squares and circles are data of the nu- 
merical transfer matrix. £ 2 an d & are amplified 20 times to 
increase visibility, but fi > £2 always holds. The lower left 
insert is a zoom in the two-channel region. The upper right in- 
sert shows the ratios ri,a of coupled to uncoupled localization 
lengths. The larger localization length is very significantly 
suppressed due to the coupling to a slow chain. Note the 
sharp recovery of the larger localization length beyond the 
band-edge E% ■ The analytical results coincide quantitatively 
with the numerical data anywhere except for specific anoma- 
lous energies: In the uncoupled case, E — corresponds to the 
commensurate wave vectors 4fei(2) = 2-7T. In the coupled case, 
E = and E w 0.03, 0.1 (very weak) and 0.17 correspond to 
2(fci + k 2 ) = 2tt, 3fci + k 2 = 2tt, 4fei = 2tt and 3fc 2 - ki = 2ir. 



small values of t, cf. Fig. 8). In this simulation, the re- 
orthogonalization method 27 was used, and length of the 
ladder is L — 10 7 with averaging over 10 3 realizations of 
disorder. The hopping integral in the fast chain t\ = 1 
was taken as the energy unit, and for simplicity, the two 
legs were taken to be equally disordered. One can see 
that as the coupling t increases the quantities rx,2 evolve 
and at t ^> W 2 approach the limits given by Eq. (101). 

The insensitivity of the localization lengths to weak 
couplings t -c 5E reflects the fact that the level spacing 
in the chains is bigger than the coupling between the 
chains, and thus wavefunctions typically do not hybridize 
much between the two legs. 

Moreover, the two families of curves for different dis- 
order strengths seem to collapse into two universal func- 
tions ri^{t/W 2 ), cf. Fig. 8. This scaling shows that at 
weak disorder W 1 and under resonance conditions 
E = Er, the numerical results approach the analytical 
ones already at a very small coupling t > W 2 . 

We can rationalize the scaling by defining a regularized 
mixing angle 7 instead of the bare 7 defined by Eq. (12). 
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FIG. 7: Localization lengths for the decoupled (t = 0, black) 
and coupled (t = 0.2, red) chains with identical hopping but 
substantially different disorder (Wi = 0.1 and Wi = 1.0) as 
a function of energy (detuning Se — 0). The solid curves are 
analytical results. The squares and circles are data from the 
numerical transfer matrix. The values of £1, £j an d £2 are 
amplified 20 times to increase their visibility. Without cou- 
pling Ci°VC2°' ~ 10 2 - ^ n the presence of coupling £1 is sub- 
stantially reduced, while £2 remains of the same order as its 
decoupled value. The insert shows the ratios n,2 of coupled 
to uncoupled localization lengths. Since t\ = t% there is res- 
onance at all energies, and thus dominance of the slow chain 
is expected. Note the sharp recovery of the larger localization 
length beyond the band-edge E} ■ There are visible anomalies 
at E — in both the uncoupled and the coupled case, which 
correspond to the commensurate condition &ku2) = 2n and 
2(fci + fe) = 2-7T. In the coupled case further anomalies exist 
at the energies corresponding to 3fci + k 2 = 2n, 4fci = 2ir and 
3/c2 — ki = 2n. However, they are every close to E = and 
too weak to be observed. 



From Eqs. (12) and (85), we find that 



tan 7 oc 



K(t U t 2 )(E - Er)^ 



(102) 



where n(t±,t2) is defined in Eq. (86). A natural way of 
regularizing the above result at resonance conditions is 
to introduce the disorder-induced "width" 8E oc W 2 in 
the form: 



tan 2 7 oc 



[K{t 1 ,t 2 ){E-E R )Y + 6& i 



(103) 



where SE scales as in Eq. (63). 

Using this regularized mixing angle, the resonant 
regime can be described more precisely by the condition 



t > max{K(ii,i 2 )|£ - E R \,SE}, 



(104) 



or, equivalently, 7 ~ tt/2. The observed scaling collapse 
in Fig. 8(a) suggests that in the weak coupling one might 
capture the behavior of localization lengths by replacing 
7 by 7 in Eq. (60). This indeed works, as confirmed by 
Fig. 8(c) where we replot the numerical data of Fig. 8(a) 
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(c) u = 0.2 
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il w 2 

FIG. 8: The ratio n,2 = f 1,2 /Ci 2 of coupled and decoupled 
localization lengths, obtained numerically as a function of the 
coupling constant t at resonant conditions, E = and 8e — 0. 
The two legs are equally disordered with a random potential 
box distributed in [-W/2, W/2]. (a) A slow t 2 = 0.2 and a fast 
ti = 1 leg: The smaller localization length increases slightly 
while the larger localization length decreases drastically, be- 
ing driven down by the slow leg, as the coupling constant t 
increases. The insert shows the dependence of n,a on the cou- 
pling constant t at different disorder strengths; all the curves 
collapse to a universal dependence on t/W 2 . The dashed lines 
in the insert show the analytic result given by Eq. (98), which 
is valid under the assumption t S> W 2 , Eq. (62). (b) Almost 
identical legs t2 ~ 0.8, t\ = 1: In this case the localization 
length of the slow leg marginally decreases while that of the 
fast leg marginally increases, (c) Results obtained analyti- 
cally upon replacing the mixing angle with a renormalized 
value, 7—^7. The parameters are the same as in (a) but with 
fewer realizations of disorder, and SE/W 2 w 0.3 in Eq. (103) 
was optimized by fitting to the numerical data in (a). The 
scaling collapse works very well in the weak coupling limit. 
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together with the analytical expressions, where 7 replaces 
7, and the number 5E/W 2 m 0.3 was optimized to yield 
the best fit. 

Note that the resonance condition can be broken either 
by detuning \E — Er\ 3> t or by increasing the disorder 
SE 3> t. Our analytic approach is based on the weak- 
disorder expansion and is therefore valid only in the first 
regime. 

4- Anomalies 

One can notice that all our numerical curves for ^12 
exhibit anomalies which are not predicted by the analyt- 
ical curves: small "peaks" appear at certain energies on 
both £1 and £2- These anomalies of localization lengths 
are due to the commensurability discussed in Sec. HID. 
This is not captured by the extended DMPK equation 
(68). However, we can identify these anomalous ener- 
gies with commensurate combinations of wave vectors in 
Eq. (57) (see the caption in Figs. 6 and 7). The anomalies 
for two chains with identical hopping but different disor- 
der have been observed numerically in Ref. 23. In this 
case there are three anomalous energies E = 0, t/2 and t, 
which correspond to commensurate combinations of wave 
vectors 2{k\ + k 2 ) = 2tt, 3fc x + fc 2 = 2ir and Ak\ = 2ir. 

B. Solutions at E 7^ 0: off-resonant regime 

Without loss of generality the off-resonant regime can 
be considered at Se — (for which the resonance is at 
Er = 0). A non-zero detuning Se merely drives En away 
from zero and induces an asymmetry of the r p as a func- 
tion of E — En. However, the mechanism of the crossover 
from resonance to off-resonance is qualitatively the same 
as in the case 5e = 0. 

Our analytical results for n,2 are presented in Fig. 9 
as functions of the dimensionless detuning E/t from res- 
onance. 

(i) Small detuning, \E\ <C t/K(t\,t-z) <C 1: The res- 
onance conditions are still fulfilled and the localization 
lengths are close to their corresponding values at E = 0. 
The leading order expansion around 7 = n/2 predicts 
that the ratios of localization lengths, n,2 only depend 
on E/t, but not on t/ti, 

|r„-r p (£7 = Q)|« (j^ , (105) 

as confirmed numerically in Fig. 9. 

(ii) Very large detuning, \E\ S> t/n(ti,t2) 1, r 1 ( 2 ) 
approaches 1 from below (above) like 

|r p -l|oci 2 . (106) 

When t is small this result is obtained from the leading 
order expansion of r p around 7 = or tt. 

(iii) For chains with equal hopping, t\ — resonance 
occurs at any energy and r p = 2C P is independent of E/t. 
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FIG. 9: Analytical results for n,a = £i,2/£i 2 as functions of 
E and E/t obtained from the extended DMPK equation in 
the weak disorder case 8E <C t for detuning Se = 0. The 
resonance energy corresponds to Er = 0. (a) ti = 0.2, t\ — 1. 
(b) t2 = 0.8, ti = 1. Close to resonance the r p only depend 
on the ratio E/t. 



C. Band-edge behavior 

Another interesting question to ask is what happens 
to the localization lengths around the band-edge E± or 
E% ? (See Fig. 1(a)) Especially, what is the behavior of 
the localization length of the fast-leg once we turn on 
the coupling t? The results from the numerical transfer 
matrix simulation and of the solution (75) of the extended 
DMPK equation, Eq. (68), are compared in Fig. 10. 
Two remarkable features can be observed in Fig. 10: 
(i) Near the band edge E = E£ where the system 
switches from one to two propagating channels, the larger 
localization length £1 (red curves) behaves in a singular 
way, as obtained from Eq. (75). As the energy tends to 
the band-edge E£ from below, £1 decreases to zero and 
shows a jump to a finite value for E > E% , where only 
one propagating channel exists. The numerical simula- 
tion (black circles) reproduces the same behavior, while 
the sharp recovering at E = E 2 is smeared by the finite 
disorder. This behavior is another drastic example of the 
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FIG. 10: Localization lengths as a function of energy near 
the edge of one of the bands. Here W = 0.2, t\ = 1, 
<2 = 0.5, t — 0.2 and Se — 0. The red curve is the result 
of the extended DMPK equation, while the black circles are 
data obtained by the numerical transfer matrix method. The 
quantitative agreement is significant except for four anoma- 
lous energies, e.g. E = and E « 0.1, 0.2 and 0.9. The 
correponding commesurate combinations of wave vectors are 
2(fci + k 2 ) = 2tt, 3fci + k 2 = 2tv, 4fci = 2tt and 3fc 2 - fei = 2vr. 
Near the termination of the lower band at E < E% the larger 
localization length is dramatically decreased, driven down by 
the slow terminating channel. At E > E~£ the localization 
length sharply recovers. In numerical simulation the sharp- 
ness is smeared by the finite disorder. 



dominant effect of the slow channel. It can be understood 
from the behavior of the Born cross-sections Eq. (60) . As 
we approach the band-edge from below, the rapidities of 
the two channels satisfy »i > 1)2. As a consequence, 
the cross-sections obey the hierarchy V 2 ^> V3 V\- 
Therefore, from Eq. (80) and (75), one can see that £1 is 
dominated by the largest cross section V 2 and shows qual- 
itatively the same behavior as £2- We emphasize that the 
mechanism of this suppression is different from that in 
the resonant regime. In the latter the suppression is due 
to 7 r~j 7r/2, which mixes the two effective variances % 2 
equally, while near the band-edge the suppression is due 
to the vanishing rapidity, which appears in the denomi- 
nators of the cross-sections. 

(ii) Anomlies are clearly seen in the numerical data for 
£1 at energies E = and E « 0.1, 0.2 and 0.9. The 
corresponding commensurate combinations of wave vec- 
tors are 2{k\ + k 2 ) = 2ir, 3fci + k 2 = 2n, 4fci = 27r and 
3k 2 — k\ — 2tt. 



VI. ONE-CHANNEL REGIME 

So far we have discussed the localization lengths in the 
two-channel regime, where the extended DMPK equa- 
tion (68) applies. In the one-channel regime (see Fig. 1) 
the second channel does not vanish but supports evanes- 
cent modes. In the presence of disorder a particle in 
propagating modes can be scattered elastically into these 



evanescent modes by local impurities. Thus the evanes- 
cent channel is coupled to the propagating channel by 
random potentials and may influence the transport prop- 
erties of the system. However, the effect of evanescent 
modes in the transport properties of ID disordered sys- 
tems is scarcely studied. 

Bagwell 28 studied in detail the transmission and re- 
flection coefficients in a multi-channel wire with a sin- 
gle (5-function impurity. The evanescent modes renormal- 
ize the matrix elements of the impurity potential in the 
propagating channels. The transmission and reflection 
coefficients of the propagating channels can be strongly 
enhanced or suppressed, nevertheless, depending on the 
strength of the impurity. 

The model in Ref. 28 was non-disordered but quite rel- 
evant to disordered systems. It is reasonable to argue 
that in ID disordered systems the effective disorder in 
the propagating channels is renormalized by evanescent 
modes, while the renormalization effect depends upon the 
strength of disorder. 

In the present two-leg Anderson model we specifically 
analyze the renormalization effect of the evanescent chan- 
nel in the weak disorder limit, which stands on an equal 
footing with the analysis in the two-channel case. Ac- 
tually, the special case t\ = t 2 and u\ — a\ has been 
studied analytically early on in Ref. 29. It was claimed 
that in the weak disorder limit the effective disorder in 
the propagating channel is significantly suppressed by 
the evanescent mode. As a consequence, the localiza- 
tion length defined through the transmission coefficient 
of the propagating channel is enhanced by a factor ~ 2 
compared to the value obtained if the evanescent mode is 
absent. However, this conclusion was unreliable because 
the average of the logarithm of transmission eigenvalue 
was not computed correctly. In contrast, we will prove 
that the evanescent channel is decoupled from the prop- 
agating channel to the lowest order in the effective dis- 
order xt defined in Eq. (48). The coupling between the 
two channels becomes relevant only at order xt- 



A. Transfer matrix of an elementary slice 

Without loss of generality, we assume that the channel 
r = 1 is propagating and r = 2 is evanescent (the upper 
branch in Fig. 1). A similar analysis applies to the op- 
posite choice (the lower branch in Fig. 1). Note first of 
all that a direct application of the Fokker-Planck equa- 
tion approach to the transfer matrix given in Eqs. (27) 
and (30) would be incorrect. The reason is the following: 
The weak disorder expansion of the parameters A, which 
leads to Eq. (51), is ill-defined in the one-channel regime. 
Note that the amplitude of the evanescent basis ip 2 (x) 
(see Eq. (26)) grows exponentially ~ e K 2M. Likewise, 
the elements of Sm x (see Eq. (27)) with evanescent chan- 
nel indices also grow exponentially with factors e 2K2 ' x ' or 
e 4K2 ' x '. Therefore, H^m^/eH is unbounded in the domain 
of the coordinate x, and the formal expansion of the pa- 
rameters A in disorder strength is divergent with respect 
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to the length L. 30 

In order to perform a weak disorder analysis, the basis 
of the evanescent channel should be chosen as 



V>± (x) = e TK2 7V2sinhK2, « 2 > 0, 



(107) 



which replaces the current-conserving basis Eq. (26), and 
the basis of the propagating channel is the same as 
Eq. (24) even though with t = 1. In this newly defined 
basis, the transfer matrix of elementary slice takes the 
form: (see App. D) 



m x = m + 5m x 



(108) 



with 



m = diag(l, 1, e- K2 ,e K2 ), Sm x = ^JJjj 



Sm\ 8m\ 



Sm 2 



whose blocks are 



5m\ = i 



8m\ =i 



en 



2sinfci Ve 



£12 



_2 _ e -i2kix 
i2k\x ^ 



2-\/sinfci sinh n 2 

g ^2 g ik\x g^2g ik\x 



5m\ 



e -K,2 e ik!X 



£21 



2-v/sin/ci sinh k 2 \ e 



^ikix 
ik\x 



— ik\x 
-ikix 



r 2 e 22 

= „ , , 

2sinhK 2 ^ £ 



-K2 _p«2 



(109) 



(110a) 



(110b) 



(110c) 



(HOd) 



where m and 5m x are the disorder-free and disordered 
part of the elementary slice m x . The transfer matrix of 
a bulk with length L is still defined by the products in 
Eq. (30). Two important points should be emphasized: 

(i) Compared with Eq. (27) in the two-channel case, 
the second and third rows and columns of Eq. (108) 
have been simultaneously permuted. The diagonal blocks 
Sm\ and Sm 2 represent the scattering in the propagating 
and evanescent channel, respectively and the off-diagonal 

1(2) 

blocks <*>w 2 (i) re P resen t the scattering between the two 
channels. In each block, the first and second diagonal 
element describe the scattering inside right- (+) and left- 
(— ) branch respectively and the off-diagonal elements 
describe the scattering between the two branches. For 
instance, 5m\+ labels the 21-element of 5m\ and stands 
for a scattering event from the left evanescent channel to 
the right propagating channel. 

(ii) The disordered part Sm x does not contain exponen- 
tially growing and/or decaying terms, and hence H^m^/eH 
is uniformly bounded for any x. Instead, the disorder- 
free part m, which is still diagonal but not unity any 
more, contains the growing and decaying factor of the 
evanescent mode per lattice spacing. The exponentially 
growing and decaying characteristics of evanescent modes 
are represented in the products of the disorder-free part 
IIs=im. 



B. Weak disorder analysis of Lyapunov exponents 

In order to calculate the transmission coefficient of 
the propagating channel, through which the localization 
length is defined (see Sec. VIC), we have to know the 
Lyapunov exponents of M(L) in Eq. (30). We are go- 
ing to determine the Lyapunov exponents by the method 
introduced in Rcf. 31. 

The Lyapunov exponents of the present model can be 
computed via the following recursive relations for the four 
vectors V^=i ... .4: 



x+l 



mxVi,x, (111a) 



Vi.x+i — rxv x V itX - 



Vj, x +i ■ (mxVi,x) 

Vj,x+1 ■ Vj,x+1 



2<i<4. (111b) 



Note that the vectors are orthogonalized by Gram- 
Schmidt procedure after every multiplication by the 
transfer matrices (108). The Lyapunov exponents are 
extracted form the growing rate of the amplitudes of the 
respective vectors: 



7i = lim — ( In L \ , 
l^2L\ \ Vhl \ 2 r 



1< i < 4, 



(112) 



in which (•) is the average over realizations of disorder 
along the strip. Moreover, {7^} are in descendant order: 



7i > 72 > 73 > 74- 



(113) 



The initial vectors V^i of the recursive relations (111) 
can be randomly chosen but must be linearly indepen- 
dent. In the absence of specific symmetry constraints the 
Lyapunov exponents are non-degenerate in the presence 
of disordered part of m^. Additionally, because of the 
symplecticity of rh x represented in Eq. (21) the Lyapunov 
exponents are related by 



73 



-72, 74 



-7i, 



(114) 



which is proved in App. D. Therefore, only the first two 
recursions in Eq. (Ill) are needed. 

In the absence of disorder the four Lyapunov exponents 
take the values: 

7l|e=0 = K2, 72| e =0 = 73|e=0 = 0, j4\e=0 = ~K 2 , 

(115) 

in which the two Lyapunov exponents corresponding to 
the propagating channel are degenerate. Therefore, we 
make an ansatz on the first two vectors, which separates 
their "moduli" and "directions", 




V 2 ,x = v 2 - 



/p(x) 
q{x) 
h{x) 

\ti{x)j 



(116) 
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in which 

si,2,3(a0,t 3 ,4(a0~O(e), (117) 
\si,2,3(x)/e\ and |*3,4(a;) / e| are bounded for all x, and 



\p(x)\ 2 + \q(x)\ 2 = l. 



(118) 



Eventually, the Lyapunov exponents are determined by 
the growth rate of {vi. x }, which is easy to be realized from 
Eqs. (112) and (116). The initial vectors of Eq. (116) are 
chosen as the eigenvectors of the disorder-free part of the 
transfer matrix m (see Eq. (109)): 






W 



2,1 




(119) 



with some p(l) and q(l) satisfying |p(l)| 2 + |<z(l)| 2 = 1- 

Note that the ansatz (116) is reasonable in the sense of 
a perturbative analysis. Consider the final vectors after L 
iterations of Eq. (Ill) with the initial condition Eq. (119). 
In the absence of disorder, it is easy to obtain V\x = 
e K 2 Ly i i anc | y 2 L — y 2 1 _ Q n t p f if weak enough 

disorder will induce perturbative effects: the direction of 
V\x will deviate from V\_\ perturbatively in the strength 
of disorder. This is characterized by the smallness of 
si,2,3(L). In other words, the exponential growth of |Vi,l| 
is dominated by m. Simultaneously, the degeneracy of 
the second and third exponents are lifted perturbatively 
As a consequence, 72 > and V2,l become exponentially 
large because of the constraint in Eq. (114). p(L) and 
q(L) are in general very different from their initial values 
p{l) and q(l), while t 3 ^{L) will be shown to remain small 
quantities of order e. 

The orthogonality between Vi tX in Eq. (116) gives 

U(x) + sl(x)p(x) + s* 2 (x)q(x) + s* 3 (x)t s (x) = 0, (120) 

in which the first three terms ~ 0(e) and the last term 
<~ 0(e 2 ). Up to the first order in disorder strength, the 
recursion (111a) gives 



«i,s+i = [e K2 + Sm 2 2 Z + 0(e 2 )] , 



(121a) 



«i,x+i»i(a: + 1) = v ltX [ Sl (x) + 5ml± + 0(e 2 )] , (121b) 
vi, x+ is 2 (x + 1) = t>i, x [s 2 (x) + 5m\z + 0(e 2 )] , (121c) 



vi, x+ \s 3 {x + 1) = [e K2 s 3 (x) + Sml + + 0(e 2 



The recursion (111b) gives 

fp(x + 1)\ 
V2 ^ (q(x + l)) =V2 '* 



(121d) 



(122a) 



V2,x+lh{x + 1) 

= v 2 , x [e- K H 3 (x) + 5m 2 +p{x) + 5m 2 tq{x) + 0{e 2 )] , 

(122b) 



V2,x+lti{x + 1) 

= -V 2 ,x [sl(x + l)p(x) + s* 2 {x + l)q(x) + 0{e 2 )] . 

(122c) 

It can be verified that the higher order terms <~ 0(e 2 ) do 
not involve exponentially growing factors, which is guar- 
anteed by the Gram-Schmidt re-orthogonalization proce- 
dure in the recursive relations (111). 

We draw two important observations from Eqs. (121) 
and (122): 

(i) The ansatz (116) is consistent with the perturba- 
tive expansion of the recursions (111). Here the con- 
sistency means that |si, 2 ,3(a0/e| and |*3,4(a;) / e| arc uni- 
formly bounded after any number of iterations, and the 
first two Lyapunov exponents can be extracted from Vj tX . 

(ii) Up to linear order in disorder strength, the recur- 
sion (121a), which determines the first Lyapunov expo- 
nent 71, is decoupled from the recursion relation (122a), 
which determines the second Lyapunov exponent 72. 
However, the coupling terms are present in higher or- 
der terms. This implies that to the leading order effect in 
disorder the evanescent and propagating channels evolve 
independently, the entanglement between the two chan- 
nels being a higher order effect. 

From Eq. (121a) one can easily calculate the first Lya- 
punov exponent to linear order in the effective variances 



71= lim - (lnTT |e K2 + Sm 2 2 Z(x)\) 
L^ooL\ l = li " '/ 



K2 



In 



1 + 



C22 



2 sinh k 2 



1 ( 2 • 

K 2 - . - 2 — ( Xi sin" - + X ~2 cos 
8 sinn k 2 V 



^4 1 

2 



2 _4 7 

2 



(123) 



in which 7 is the mixing angle defined in Eq. (12). The 
minus sign of the leading order corrections implies that 
the first Lyapunov exponent is reduced in the presence of 
weak disorder. 

Eq. (122a) is exactly the same as in a single chain An- 
derson model, for which the Lyapunov exponents are al- 
ready known. 5 The second Lyapunov exponent takes the 
value 



l2 ~2V 1 + 0(xi), 



(124) 



where V\ is the Born cross-section given in Eq. (60). 

Eqs. (123) and (124) are our main results for the one- 
channel case, yielding the localization length and the 
renormalized decay rate of evanescent waves. 



C. Localization length and evanescent decay rate 

The two Lyapunov exponents calculated above can be 
identified in transport experiments. In general a two- 
probe experiment has the geometry of the form "lead- 
samplc-lead" , in which the two leads are semi-infinite. 
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The current amplitudes (not the wave amplitudes) are 
measured in leads. In the propagating channels both 
right (+) or left (— ) modes exist in both of the leads. 
However, the situation is rather different in the evanes- 
cent channels: There are only growing modes (— ) in the 
left lead, and only decaying modes (+) in the right lead. 
These modes do not carry current at all. 28 ' 32 Hence the 
current transmission and reflection coefficients are only 
defined in propagating channels regardless of the wave 
amplitudes in evanescent channels. In terms of the trans- 
fer matrix M(L), this restriction on the evanescent chan- 
nel implies that 



fat(L)\ 

a i ( L ) 
a+(L) 

V o J 



Ml Ml 
Ml M| 



/Mi)\ 

Mi) 
o 

Way 



(125) 



From the scattering configuration (125) one can de- 
rive an effective transfer matrix for the propagating chan- 
nel. The evanescent amplitude a^~(l) can be expressed in 
terms of the propagating amplitudes as 



%(!) = - 



Mi 



[M^a+CO+M^Ml)]. (126) 



Substituting Eq. (126) into Eq. (125) we obtain 



a+(L) 
ai(L) 



ML) 



Mi) 



in which the elements of X(L) take the form: 



x+ = mix + AM lii AM U 



X+ = Ml+ + AMl+, AM 1 



i+ 



Ml+Ml~ 
Ml~_ 



M\+Ml~_ 



Mi 



x- 



x- 



x- 



X 



+* 



(127) 



(128a) 



(128b) 



(128c) 



X(L) is the effective transfer matrix for the propagating 
channel. Note that its elements are modified form the 
values in the absence of the evanescent channel. One can 
easily verify that X(L) satisfies time-reversal invariance 
and current conservation conditions as (29) in the single 
chain case: 11 ' 12 



X* = (7lX(Tl, XV3X 



03 • 



(129) 



However, X(L) does not evolve multiplicatively with the 
length L any more. The transmission coefficient is deter- 
mined through X(L) in the usual way 11 ' 12 



T(L) = \X+\- 



(130) 



where X+ is defined in Eq. (128a). 

Eqs. (128) and (130) exactly determine the transmis- 
sion coefficient of the propagating channel. A full so- 
lution requires extensive calculations. However, if the 



disorder strength is weak, as analyzed in Sec. VI B, the 
coupling between the two channels is small, so that the 
contribution of the evanecent channel, AM 1 ^" is negli- 
gible. Indeed, from Eqs. (Ill) and (116), using initial 

vectors V M = (0 l) T and V 2 ,i = (l 0) T , re- 
spectively, we can extract the various matrix elements of 
M(L), in particular 



AMl+(L) _ Sl (L)t 4 (L) 



Mlt(L) 



P(L) 



0(e 2 ). 



(131) 



This proves that the contribution of the evanescent chan- 
nel is subleading at weak disorder. To leading order the 
transmission coefficient is simply given by the propagat- 
ing channel as 



T(L) 



\Mlt\ 



\v 2 , L p(L)\- 



From this the localization length is obtained, 



lim ±(\nT(L))c2V 1+ 0(xt) 



8 sin 



in fci V 2 2 / 



(132) 



(133) 



Eq. (133) implies that to leading order in \1 the localiza- 
tion length in the propagating channel equals the inverse 
of the second Lyapunov exponent obtained in Eq. (124). 

Similarly to Eq. (91), we can introduce the localization 
length enhancement factor 



(0) 



(134) 



in which £^ is the localization length of the leg 1 (with 
the larger hopping) in the absence of inter-chain coupling. 

On the other hand, the inverse of the first Lyapunov ex- 
ponent in Eq. (123) should be associated with the evanes- 
cent decay rate which is slightly modified by disorder. 

The analytical results (133) and/or (134) are compared 
with numerics in Figs. 6, 10 and 11. Figs. 6 and 10 corre- 
spond to the weak coupling case t < t c (see Fig. 1(a)) and 
Fig. 11 to the strong coupling case t > t c (sec Fig. 1(b)). 
The remarkable agreement confirms the weak disorder 
analysis developed in this section. 

We specifically analyze the typical behavior of the en- 
hancement factor r\ (E) close to the band-edge in the 
case of t < t c , where the system switch from one to two 
propagating channels. From Eqs. (12) and (18) it is not 
hard to obtain: at the band edge E^ , when coupling is 
weak n(-E^) deviates from 1 like 



i-n(E+)(x 



El 



(135) 



If E is away from E£ , i"i(E) increases linearly, i.e., 

r 1 {E)-r 1 {E+) oct 2 (E-E+), (136) 

with a fixed but weak coupling t. A typical curve for 
r\(E) is shown in the upper right insert in Fig. 6. 
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FIG. 11: Localization length as a function of energy in the 
one-channel case. Here t\ = 1, tz = 0.1, 5e = and the 
amplitude of disorder is W = 0.1. The solid curves are ana- 
lytical results. Black curves correspond to uncoupled chains. 
The red one corresponds to the upper polariton (conduction) 
band (propagating channel) for strong coupling t = 1, which 
is obtained by omitting the lower polariton (valence) band 
(evanescent channel). The squares and circles are data of 
the numerical transfer matrix. The quantitative agreement is 
significant except for the anomalous energy E ~ 1.0, which 
correspands to the commesurate combination of wave vectors 
4fci = 2-7T. The coincidence between analytics and numer- 
ics confirms that the evanescent channel is decoupled to the 
propagating channel in weak disorder limit. 



VII. SHAPE AND POLARIZATION OF THE 
WAVEFUNCTIONS 

In certain applications, such as exciton-polaritons, the 
two linearly coupled types of excitations (represented by 
the two chains) are very different in nature. This makes 
it in principle possible to probe the original excitations 
separately from each other. For a two-leg atomic chain 
one can imagine probing the amplitude of wave function 
on each of the spatially separated legs. For polaritons the 
analogue would be a separate probing of cavity photons 
or excitons, e.g. by studying the 3D light emitted due to 
diffraction of cavity photons at surface roughnesses or by 
studying the exciton annihilation radiation or the elec- 
tric current of exciton decomposition provoked locally. 
Therefore it is of practical interest to be able to manip- 
ulate the strength of localization of one of the original 
excitations by coupling them to the other. 

A. Numerical analysis 

With this goal in mind we have carried out a numerical 
study of the amplitude of wave functions on either of the 
chains in each of the distinct parameter regimes discussed 
above. We numerically diagonalized the Hamiltonian (1), 
cf. Fig. 12, choosing t x = 1, t 2 = 0.2, W = 0.4, t = 0.04. 
The length of the ladder was taken to be L = 10 3 and 



periodic boundary conditions were used. With the above 
parameters the localization lengths of the decoupled legs 

were of the order of ^ ~ 10 3 and £2 ~ 10' f° r energies 
close to the band-center. In Fig. 12 the black curves 
depict the amplitudes of eigenfunctions on the fast leg 1, 
while the red curves show the corresponding amplitudes 
on the slow leg 2. 

Our main findings are: 

(i) E~2 < E < E^: The energy is far from resonance, 
and only one channel exists. As shown in Fig. 12(a), most 
of the weight is on the fast leg. The amplitude on the slow 
leg is small but the spatial extension of the component ip2 
is the same as that of ip\ on the fast leg, which is almost 
unaffected by the chain coupling. Thus the coupling can 
create a nonzero amplitude on the chain 2, in the energy 
region where the decoupled chain 2 cannot support any 
excitations. The spacial extension is controlled by the 
localization properties of the leg 1. 

(ii) t/n(t\,t2) < E < E^: The energy is in the two- 
channel, off-resonant regime. The wavefunction compo- 
nents ipi and ^2 are characterized by both localization 
lengths £1 and £2- However, the relative weights of the 
parts of the wavefunction with the smaller and the larger 
localization lengths fluctuate very strongly from eigen- 
state to eigenstate. This is shown in Fig. 12(b) and (c), 
with two adjacent energy levels, which were properly se- 
lected. In Fig. 12(b), ip2 consists almost entirely of a 
component with the smaller localization length, while the 
fast leg clearly shows contributions of both components. 
In Fig. 12(c), both ipi and il>2 consist almost entirely of a 
component with the larger localization length. In brief, 
the former can be thought of as a state on leg 2, which 
weakly admixes some more delocalized states on leg 1, 
while the latter wavefunction is essentially a state of leg 
1 which admixes several more strongly localized states on 
leg 2. 

We have checked in specific cases that this interpre- 
tation is indeed consistent (see Sec. VII B): In the off- 
resonant regime the wavefunctions can be obtained per- 
turbatively in the coupling t, confirming the picture of 
one-leg wavefunctions with small admixtures of wave- 
functions on the other leg. Off resonance, the pertur- 
bation theory is controlled even for appreciable t, since 
the matrix elements that couple wavefunctions of similar 
energy are very small due to significant cancellations aris- 
ing from the mismatched oscillations of the wavefunctions 

(fci — fc 2 > £2 ) on the two legs. Resonance occurs pre- 
cisely when at a fixed energy k\ — k 2 becomes too small, 
so that the modes on both legs start to mix strongly. 
A closer analysis of the perturbation theory in special 
cases shows that the perturbative expansion is expected 
to break down at the resonant crossover determined fur- 
ther above. 

(iii) \E\ < t/ n(ti,t 2 ): If the energy is in the resonant 
regime, the two localization lengths are of the same or- 
der £1 ~ £2 ~ £2°' and the spatial extension of both 
wave function components is governed by the localization 

length of the decoupled slow chain. This is illustrated 
in Fig. 12(d). 



23 




FIG. 12: Typical wavefunctions in different regimes. Pay attention to the varying scales for the amplitudes in the various cases. 
Parameters are ti = 1, ta = 0.2, W = 0.4, t — 0.04, while E selects the regime. The length of the ladder is L — 10 3 with 
periodic boundary condition. The black(red) curves are the amplitudes on the l(2)-leg. (a) E = 1.5955, at which only one 
channel exist, (b) and (c) off-resonant regime: E = 0.3250 and E — 0.3242 are a pair of adjacent levels of "opposite type", (d) 
E ~ 0, which is within the resonant regime. 



B. Perturbative analysis 

The properties of eigenstates at different energy 
regimes can be explained by applying a perturbative anal- 
ysis on the coupling t. First, we define the relevant quan- 
tities of decoupled legs as follows: the eigenstate of the 
u-leg with eigenenergy E vn is ij) vn (x). The corresponding 
localization length is $\ where we assume £^ 3> £2°^ 
in order to reveal the resonance-off-resonance crossover. 
The mean level spacing inside the localization volume is 
A„. Because in one dimension a particle is nearly ballis- 
tic in its localization volume, which means its wavevector 
is nearly conserved and its amplitude is almost uniform, 
we introduce a simple "box" approximation on the eigen- 
states as the following: Inside the localization volume, 

1 




ik v X 



(137) 



up to a random phase, in which and k T are the lo- 
calization length and the wavevector at the energy E vn . 



Outside the localization volume ip vn (x) = 0. 

Now we turn on a weak enough coupling t and calculate 
the deviation of an energy level E\ n on the 1-leg. Up to 
second order in t, the deviation is 



SE {2) 

0Cj ln 



* 2 E 



|J dxipl n (x)^ 2m (x)\ 



E 



2,1, 



(138) 



(2) 

In order to estimate the value of SE{^ by the r.h.s. of 
Eq. (138) we have to make clear three points: 

(i) The summation is dominated by the terms with the 
smallest denominators, whose typical value is the mean 
level spacing A2. 

(ii) The typical value of the integral on the numerator 
can be estimated by the "box" approximation introduced 
above, which gives 



dxipl n ip 2m 



« 0) 



(139) 
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(iii) We should consider more carefully how many dom- 
inant terms there are in the summation. We can easily re- 
alize that a state ipi n {%) on the 1-leg can couple to about 

Ci ''/^ '' states ip 2m (x) on the 2-leg. However, the value 
of the summation is different from a naive deterministic 
evaluation because the random signs of the denomina- 
tors. If we neglect the correlation of these random signs, 
according to the central limit theorem, the fluctuation of 

SB® * 



5E\ 



(2) 



ft 



(0) 



(140) 



The validity of the perturbation analysis is guaranteed if 



SE\ 



(2) 



< Ai 



(141) 



which means there is no level-crossing in the localization 
volume of the 1-leg. To estimate the relevant quantities in 
Eq. (141), for simplicity we assume t\ ^> t 2 and a\ = erf. 
If the energy E = E ln is close to the resonant energy En, 
according to Eq. (10) and (14), we obtain 

|*i - k 2 \ ~ \E - E R \ (h - t 2 )/tit 2 . (142) 
The mean level spacings are 

Av~t v /£l°\ (143) 
and the localization lengths satisfy 

~ tl/tl (144) 

Substituting Eqs. (142), (143) and (144) to Eq. (141) we 
obtain the condition 



t< \E-E R \(t 1 -t 2 )/h, 



(145) 



which is consistent with the result of Eq. (85) with 
t\ ^> t 2 . Therefore, Eq. (141) is essentially equivalent 
to the criterion for being off-resonant (AV > 0) at weak 
coupling t. 



VIII. LIMIT OF VANISHING HOPPING ON 
THE "SLOW" LEG 

In the present work we are particularly interested in 
the case where the localization lengths of the uncoupled 
legs are parametrically different ^ 3> £2°^- Accordingly 
we refer to the two legs as the "fast" and the "slow" 
one, respectively. So far we have analyzed the model 
extensively in the limit where the disorder is weak on 
both legs and thus > 1. 

Another interesting situation is the case where the hop- 
ping strength on the slow leg vanishes t 2 = 0, or is weak 
enough. This is experimentally relevant for polariton sys- 
tems in which the exciton hopping is weak as compared 
to the disorder potential. In this case the dimension- 
less disorder parameter which we introduced previously 



diverges x\ ~~ 00 > an d formally £. 



(0) 



even if the dis- 



order strength on this second leg is arbitrarily small. For 
this reason the perturbative analysis in both of xl breaks 
down. Nevertheless, this limiting case can be solved ex- 
actly, too, but requires a different treatment which goes 
beyond the previous weak disorder analysis. 

If t 2 = 0, the second leg is composed of mutually non- 
connected sites, which form a comb structure together 
with the first leg. The Schrodinger equation (4) takes 
the form 
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+ 1) + V(x 

E - e x i 

t E- 
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where 



tf(a;) 



V>i(z) 

1p2{x) 



(146) 



(147) 



describes the amplitudes on the two legs, respectively. 
We obtain the effective Schrodinger equation for ip\(x) 
by eliminating ip 2 (x) in Eq. (146): 



-tl [lp!(x + 1) + 1pt(x 

where 



1)] = (£-e xl )Vi(z-l), (148) 



E - Se - e x2 



(149) 



Note that e x \ has the meaning of an effective disorder po- 
tential on leg 1. Furthermore, if \e x2 \ <C \E— 5e\ Eq. (149) 
can be expanded as 



Cxi 



E-5e 



+ 



t/ 



(E - 5e)- 



;£x2 



+ 0(e 2 ). (150) 



The first term on the r.h.s of Eq. (150) is a homogeneous 
potential shift. The second term is an effective disorder 
potential of zero mean. 

Eqs. (148) and (150) represent a single-chain problem, 
which can be solved exactly. The dispersion relation of 
the disorder-free part is determined by 



-2ti cos A; + 



E, 



E-Se 

which gives the two non-overlapping bands 



(151) 



Se 

E T {k) = -tiCosfc+ — - (-1) 7 



ti cos k + 



Se 



+ t 2 . 



(152) 

Of course, this coincides with Eq. (15) for t 2 = 0. Using 
the result for a single chain Anderson model 5 we obtain 
the localization length as 



xi 



xt tan 2 7 



i sm 
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+ 0{xi,xi)- 



(153) 



Here, the disorder on the leg 2 is measured by the dimcn- 
sionlcss ratio 



X2 



4 



(154) 
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FIG. 13: Localization length as a function of energy for three 
values of £2- £i = t = 1 and W\ = W2 — 0.1 are kept fixed. 
t2 = 0.2, 0.05 and 0.001 capture the weak disorder limit, the 
intermediate regime and the limit of a nearly disconnected 
slow leg, respectively. The symbols are data of a numerical 
transfer matrix calculation. The solid and the dashed curves 
are the analytical results from Eqs. (133,153). The agreement 
with numerics in the limiting cases is very good. The localiza- 
tion length increases monotonically with for a fixed energy. 
An anomaly due to the commensurate wave vector 4k = 2tv 
appears at E w 1.0. 



7 is the mixing angle denned in Eq. (12) with £2 = 0. 
Comparing the result (153) with Eq. (133), which de- 
scribes the one-channel case in the weak disorder limit 
(W2 <C £2), one sees that the two limits do not com- 
mute. This is similar to the non-commutativity of the 
limits of weak disorder and weak inter-chain coupling in 
the resonant case. Similarly as Eq. (63) the characteris- 
tic disorder energy scale is the mean level spacing in the 
localization volume when £2 = 0, which can be estimated 
as 



§E ~tx/€~ ma x{oi,al}/ti. (155) 



A perturbative analysis in £ 2 is valid only if £2 <C SE. We 
expect a crossover to the regime of strong hopping on the 
slow leg when £ 2 ~ 5E. 

The non-commutativity of the two limits is illustrated 
in Fig. 13 where we compare the two analytical limits 
with numerical simulations at fixed disorder W\ = W2 = 
0.1 and hoppings t\ = £ = 1. Three values of t 2 are se- 
lected to cover the crossover from the weak disorder limit 
(£2 ^> W2) to the limit of a slow leg with disconnected 
sites (£2 — > 0). Eqs. (133) and (153) are indeed seen to 
capture the two limits very well. Note that the local- 
ization length increases monotonically with £ 2 for a fixed 
energy, as one may expect. 



IX. CONCLUSION AND POSSIBLE 
APPLICATIONS 

The most important potential application of our the- 
ory is in the realm of polaritons in quasi-one dimensional 
semiconductor structures 14,15 . Here the fast chain corre- 
sponds to the electromagnetic modes ("light") confined 
in a one-dimensional structure and therefore having a 
parabolic dispersion with a very small mass (large £1) 
at small wave vectors. The slow chain corresponds to the 
Wannier-Mott excitons, an electron-hole pair coupled by 
Coulomb attraction. The mass of the exciton is typically 
10 4 times larger than that of "light". Surface roughness 
of the one-dimensional structure and impurities therein 
produce a disorder potential acting on both excitons and 
" light" 16 . Experimentally, one can easily probe the inten- 
sity of the " light" component by measuring the intensity 
of 3D photons that emerge due to diffraction from the sur- 
face roughness. The amplitude of the wavefunction of the 
exciton's center-of-mass is more difficult to access, but in 
principle still possible, e.g. via stimulated or spontaneous 
exciton recombination and the related radiation. Another 
application is related with one-dimensional structures in 
cold-atom traps. By this technique one can construct and 
study coupled one-dimensional chains in the same way as 
it was recently done for a single chain 34 ' 35 . 

While all possible regimes can be achieved in a 
system of cold atoms, the most relevant regime for 
one-dimensional cavity polaritons is that described in 
Fig. 1(b), where for each band (the lower and the up- 
per polariton bands) only one channel exists. This is 
because the exciton-light coupling is typically stronger 
than the narrow bandwidth £2 ~ l/m cxc of the excitons. 
The localization length in the case relevant for the upper 
exciton-polaritons is shown in Fig. 11. 

As expected, the localization length tends to zero near 
the bottom and the top of the band, due to the vanish- 
ing rapidity. Note that in the context of polaritons the 
upper band edge does not exist, since the light has an 
unbounded continuous spectrum. In our model the top 
of the band appears merely due to the discreteness of the 
lattice. In the center of the band the localization length 
is of the order of the localization length for the uncoupled 
"light" component. More important is the distribution of 
amplitudes of " light" and " exciton" wavefunctions which 
are similar to Fig. 12(a), with "light" being represented 
by the wave function tpi on the fast leg and the "exciton" 
part being represented by the wave function ip2 on the 
slow leg. One can see that the coupling to light makes the 
exciton wave function spread over a distance £ ~ £^ of 
the order of the localization length of light. This is much 
larger than the maximum exciton localization length £g 
in the absence of coupling. The price for the "fast tran- 
sit" is that the amplitude of the exciton wave function is 
small. This means that the transfer of a locally created 
exciton to distances of order £1' is possible, but occurs 
with reduced probability. 

In the long search for light localization (see the paper 
by A. Lagendijk, B. van Tiggelen and D. S. Wiersma in 
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Ref. 2) the crucial point was to achieve a smaller localiza- 
tion length of light. Our results show that this can also 
be achieved by coupling light to excitons near the bottom 
of the upper polariton band. 



We would like to emphasize that the model considered 
above does not take into account an important property 
of polaritons, namely their finite lifetime due to recombi- 
nation of excitons, and the out-coupling of the light from 
the waveguide. This limits the coherence of polaritons 
and inhibits Anderson localization. More precisely, the 
effects of Anderson localization are only relevant if the 
time to diffuse up to the scale of the localization length 
(which by Thouless' argument is of the order of the in- 
verse level spacing in the localization volume) is smaller 
than the life time of the excitons. Further crucial aspects 
are interactions among polaritons at finite density, and 
the related possibility of interaction-induced derealiza- 
tion and Bose condensation of polaritons. 36,37 A complete 
theory of localization of hybrid particles like polaritons 
should take into account all these issues. 39 



Let us finally discuss the role of dimension for our re- 
sults. We have found that under resonant conditions the 
localization lengths of two coupled chains are of the order 
of the localization length of the more localized, uncoupled 
leg. We may interpret this phenomenon as a manifes- 
tation of the fact that in Id the mean free path is the 
relevant length scale that sets the localization length. In- 
sofar it is not surprising that the backscattering rate, and 
thus the "worst" leg of the chains determines the localiza- 
tion properties of a coupled system. However, the close 
relation (proportionality) between mean free path and lo- 
calization length is special for one-dimensional systems. 
In contrast in two dimensions the localization length be- 
comes parametrically larger than the mean free path at 
weak disorder. In d > 2 most eigenstates are even de- 
localized in weak disorder. Accordingly, we expect that 
the localization length is not so simply determined by 
the properties of the more disordered part among two 
coupled systems. Nevertheless, since the proliferation of 
weak-localization and backscattering leads to complete 
localization also in 2d (in the absence of special symme- 
tries), we expect that a well propagating channel becomes 
more strongly localized upon resonant coupling to a more 
disordered channel. This may apply, e.g., to 2d polari- 
ton systems. However, in higher dimensions d > 2 such 
a coupling might have a rather weak effect. We expect 
that a "fast" channel is not affected much by a more dis- 
ordered parallel channel. That such a trend exists indeed 
at high enough dimensions can be shown in the case of 
two coupled Bethe lattices, 40 which can be viewed as the 
limit of arbitrarily high dimensions. 

We leave the investigation of problems in higher di- 
mensions, and possible implications for interacting few- 
particle problems for future work. 



Appendix A: Transfer matrix of an "elementary 
slice" in the current-conserving basis Eqs. (24, 26) 

In this appendix we derive Eqs. (27) and (29). In 
Eqs. (19, 20), m x is a symplectic matrix which by def- 
inition satisfies 

Jfa x = J, (Al) 

m x = m*, (A2) 

where 

J =(-°iS)' ^ 

and rh x is the matrix obtained from m x by transposition. 
Define the new matrix m x 

m x = U^ 1 m x U x , (A4) 

where the rotation matrix is 

u - s («£, <_ ) ■ < A5 > 

with ol x defined by Eq. (28). The corresponding inverse 
matrix is given by: 

u " = 1 (-1;', «!:) < A6 > 

where A is the diagonal matrix 

A = a x a.* x _ 1 — a* x a x -\. (A7) 

The crucial point is that by current conservation Eq. (23), 
A is independent of coordinates and is proportional to the 
unit matrix in channel space: 

A = il. (A8) 

Note also that, by construction, the rotation matrix XJ X 
obeys the disorder-free Schrodinger equation Eq. (19): 

U x+ i = ih x \ ?=0 U x . (A9) 

It follows immediately from Eqs. (A9, A4) that in the ab- 
sence of disorder m x = 1. In the presence of weak disor- 
der the matrix m x acquires a small coordinate-dependent 
correction proportional to e x which is given by Eq. (27). 

Next, by inverting Eq. (A4) and plugging into Eq. (A2) 
one obtains: 

sWm.sW^m* (A10) 
Using definition of a x Eq. (28) one can readily show that 

s« = (u;)- 1 ^ = Si (Aii) 

is real and independent of x. This immediately reduces 
the time-reversal symmetry condition Eq. (A10) to the 
form in Eq. (29). 
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The same procedure applied to the symplecticity rela- 
tion Eq. (Al) results in the following constraint (using 
m£ = mt): 



(A12) 



where 



s! 3 > = Ut JU, 



-AS 3 = -iS 3 (A13) 

is independent of coordinate due to current conserva- 
tion. Thus we obtain the current conservation condition 
in Eq. (29). 



S(F + <5F)S T = F + F^ + F^ 2 \ (B7b) 

We have introduced perturbation terms on the r.h.s. of 
the two equations as 



F (1) = -FA L + A L F + S L F + FS 



Li 



(B8a) 



F( 2 ) = -A L FA L + £ L FE£+A L FE£-i: L FA L , (B8b) 



F (1) = -FA^ + A L F + S fj F + FS, 



(B8c) 



Appendix B: Perturbative calculation of SX up to 
second order 



F< 2 ) =-A L FA|, + S L FS L +A L FS L -S L FAJ,, (B8d) 



Eqs. (38), (45) and (46) fully determine the variation of 
the eigen-system of the Hermitian matrix R. It is given 
by <5A, which characterize the "perturbation" i5R. We can 
therefore use standard perturbation theory to expand SX 
into powers of disorder on the additional slice. In this 
Appendix we calculate SX up to the second order, which 
is necessary to derive the Fokker-Planck equation (51). 

We introduce some quantities which are convenient 
to present the results. Analogously to ol x , defined by 
Eq. (28), we define 

X = <**u, (Bl) 

where u is the unitary matrix in Eq. (34). Since a x 
describes the propagation in the plane-wave basis on the 
individual chains, and u is the "polarization" matrix, we 
can consider j3 x as describing clean propagation in the 
"polarized" plane-wave basis. Furthermore, analogously 
to the blocks in Eq. (27), we can define two quantities 
on the "polarized" basis, related with the forward- and 
back-scattering of the right-moving particle off the slice: 

A x = i/3 x E x (3 x , 

S x = 0l~e x f3* x , 

which arc 2x2 matrices. It is easy to realize that A x is 
anti- Hermitian and H x is symmetric. The corresponding 
left-moving quantities are complex conjugates of them. 
The perturbative series of SX are functions of elements of 
A x and S x . For simplicity of further notations, we define 



(B2) 



and 



F = VF 2 - 1, 



AF = F 1 - F 2 



(B3) 



(B4) 



In order to facilitate the perturbative calculation, we 
adopt a parametrization of R + <5R as in (38), but with 
F -> F + SF, F -)• F + SF and u -> u + Su, in which 



SF = a/(F + <5F) 2 - 1 - VF 2 - 1. (B5) 
Substituting this into Eq. (46) and (45) we obtain two 
coupled equations for SF and the 2x2 matrix S which 
captures the incremental change of the polarization basis, 

S = 1 + u^u, (B6) 



as 



S(F + <5F)S t = F + F (1) +F ( 



(B7a) 



where F^ and F^ are linear in disorder, while F^ 2 ) 
and F^ 2 ) are quadratic. Additionally, F^ and F^ 2 ) are 
Hermitian, but FW and F' 2 ) are symmetric. 

We expand SF and Su in disorder strength. From the 
latter we calculate the corresponding variations of angu- 
lar variables. Without going into the details of the cal- 
culation, we present the results up to the second order in 
disorder. 

To first order the corrections 

«5A« = (SF^,SF^\ 50^,5^,5^,5^) (B9) 



are given by 



«f)=F« ge{l,2}, 



(BlOa) 



56W = ^Re(F^), (BlOb) 



H ( ) = o ~ ^ ^ - 6( P ( ) COS0 i ( B10c ) 



(i) = _l / imF^ ^ Imfgj ' 

2 v F 2 h 



(BlOd) 



6<p 



(i) 



AF 



ImjF^e^) csc(9, (BlOe) 



where the subscripts denote the matrix elements of the 
"perturbations" F^ in Eq. (BIO). We recall that the 
"perturbations" in Eq. (BIO) are L-dependent. 
The second order corrections 

SX^ = {5F[ 2) ,5F^,56^\5^ 2 \5^\5^) (Bll) 
are more complicated. However, we recall that our aim 



is to calculate the correlators SX^pSX^ and SX { /' in 
Eq. (51). To avoid repeating calculation, we should ex- 
press the SX^'s in terms of the first order corrections 
Eq. (BIO) as far as possible. We obtain 



v(2) 
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*F fi ( 2 ) = F< 2 ) - (-1)' 



AF ' 



<?€{!> 2}, 



69^ 



AF 



Re(Fge^) + a^59^ + * sin 29 (V 1 ')' , 



^ (2) = aL 2) + b W 6^ + b^d^ + c+e ( 2 ) - a W c^ - \ sinftW* 1 V 1 ) - cos0<M 2) , 



5ip (2) = — ImfF^e 1 ^ csc 9 + a (1) 6<p w -cot 959^6<p (1 \ 



(B12a) 

(B12b) 

(B12c) 
(B12d) 
(B12e) 



in which 



"'■sK-»f) 



(B13a) 



e (2) = * 



sin 2^ 



+ *(V (1) ^ (1) sin(9cos2V>. 



(B13f) 



a (2) = 1 ( Igjg ± ImF< 2 r 

2 \ i?2 Fi 



(1)_ 1 , 



2 VF 2 



SF^ ± -^SF : 



n 



(B13b) I n practice, we first calculate all the correlators SX^SX^ 
by Eq. (BIO). At the same time we obtain the correla- 
tors relevant for the products of first order terms on the 
r.h.s. of Eq. (B12). Finally, after evaluating the disor- 

(B13c) der average of in Eq. (B13b), 5Xf\ can be obtained. 



1 FiF 2 \ 



(B13d) 



Appendix C: co-efficients of Eq. (68) 



d« = (V W cos0 + <$V (1) , 



(B13e) The co-efficients of Eq. (68) are 
I 



Cl 



Fo 



02 = 2 AF r6 ' 



(Cla) 
(Clb) 



C3 



AF 2 



[(f 2 - f 2 - 2) (r 5 + r 4 ) + 2FxF 2 (r 5 - r 4 ) 



(F 2 + F|) (1 - u 2 ) d u T 6 



AF 2 



F 1 F 2 uT 6 cos 2^, 



(Clc) 



c 4 = 



1 

AF 



- - F 9 ~ F, - \ FiF, u 



+ ^F 2 a u r 4 + $Fid„r 5 - 2 4$ (r 6 + 2ud u r 6 ) 



F x 



AF 



sin2?/>, 



(Cld) 



29 



Cll=F?T u 



ci2 = 2FiF 2 r 3 cos2V>, 



(Cle) 
(Clf) 



IF / \ 



c 2 2 — f 2 r 2 , 



023 = 2^(^2 + ^008 2^ T 5 



Fi F-i u \ 

C24 = -F 2 | ^r 3 + 2 .' .. ^ r 5 sin 2?/-, 



Fx AF 713^2 



C33 



^3 + ^2 (^1 + P 2 + cos 2V-) T 6 



(l-« 2 ), 



1 / F\r2 Fir 2 F1F2 \ 

C34 = -~r-^ — = — T 4 H = — r 5 - 4 uT 6 sm2ip, 



AF \ F x 



AF 



(Clg) 
(Clh) 
(Cli) 
(Clj) 

(Clk) 
(Cll) 
(Clm) 



c 44 



2 
+ 



1 + 



1 F, 



2 \F l 



1 + 



1 Fo 



VT 



2 + ^)r 4+ ( 2 + AV 5 



2 Vf 2 
i 

AF J 



r 2 -r 6 



+ 



r 3 + 1 + 



lf\F 2 _ju 1 / F X F 2 F X F 2 \ u 2 F X F 2 

2F 1 F 2 3 VT^&f{ F 1 1 F 2 2 l-u 2 AF 2 



F? + F 2 
AF 2 

cos 2ip, 



(Cln) 



where the F e 's are the two diagonal elements of the ma- ties introduced in Eq. (C) are defined below. Notice first 
trix (B3) and AF is defined in Eq. (B4) . The new quanti- that The T n are functions of u defined by 



ri (2 )(«) = V 1 + V 2 + W 3 - (+)2 (V 2 -V 1 )u + (Vi + V 2 - 4V 3 ) u 2 , (C2a) 

T 3 (u) = (V 1 + V 2 -W 3 )(l-u 2 ), (C2b) 

r 4 (5)(«) = [^1 - ^ 2 + (-) {Vx +V 2 - W 3 ) u] (1 - u 2 ) , (C2c) 

r 6 ( M ) = y 1 + F 2 -(yi + y 2 -4y 3 ) u 2 , (C2d) 



where V±, V2 and V 3 are the three Born cross-sections limit L » 1, we need the values of c 3 and C33 in the limit 
defined in Eq. (60). In order to solve Eq. (78) in the 
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The transfer matrix of an elementary slice (108) is 



lim c 3 = {\V 1 -V 2 \ -d u T e )(l-u 2 ), 



=U x + 1 m x U ;i 



Sm x , 



(D3) 



lim c 33 - (V3 + r 6 ) (1 - u 2 ) , 



which is Eq. (79). 



Appendix D: Transfer matrix of an "elementary 
slice" in the basis Eqs. (24, 107) 

The derivation of Eq. (108) goes in parallel with the 
derivation of Eq. (27) in App. A. However, the rotation 
Ua; is constructed in such a way that 



U, 



/ ^+{ X ) vr(x) \ 

V 2 + (i) V> 2 ~(i) 

ipt(x- 1) ipi(x-l) 

V v 2 + (o) ^wy 



(Dl) 



where ipf 2 (. x ) are defined by Eqs. (24) and (107). Com- 
pared with Eq. (A5) for the two-channel case, the second 
and third columns of Eq. (Dl) have been permuted, and 
the columns corresponding to the evanescent channel are 
coordinate- independent. The inverse of Ua; is 



(-iipi{x-l) iipi(x) 
iiP+(x-l) -ii>t( x ) 

-v 2 "(0) v 2 "(i) 



V 







tf+(0) 







-v4(i)/ 



(D2) 



where m and Sm x take the form given in Eq. (109). Fol- 
lowing the same procedure as in App. A, one can also ob- 
tain the symmetry constraints on the matrix m x , which 
is imposed by the reality and symplecticity of the matrix 
m x . Without going into details we present the following 
results: The reality relation (A2) gives 



m* = Aim x Ai, Ai = 



<7i 

1 



(D4) 



The symplecticity relation (Al) gives 



m K A 3 m x = A 3 , Ai = ^ 3 . 



(D5) 



Finally, it is not hard to show that the Lyapunov expo- 
nents of the products (30) satisfy the symmetry property 
stated in Eq. (114). 
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